pynx.ptycho
: 2D Ptychography#
Description#
This modules allows the simulation and analysis of ptychography experiments, with the following features:
2D ptychography using a variety of algorithms: alternating projections, difference map (Thibault et al.), maximum likelihood conjugate gradient
Works with any type of illumination (sharp or not)
Object and/or Probe reconstruction
Probe analysis (modes and focus)
Incoherent background optimization
GPU implementation using CUDA/OpenCL is available (and recommended), and is the main focus of current development:
example speed on single V1OO GPU card as of 2019/06: 13 ms/cycle for 1000 frames of 256x256 pixels and a simple alternating projection (17 ms/cycle for DM and 34 ms/cycle for ML)
GPU implementation allows using modes for probe and object
Maximum likelihood minimization (Poisson noise, regularisation)
simple usage scripts are provided to analyse data from CXI files, ESRF beamlines (id01, id13, id16a), and ptypy files.
API Reference#
Note that the Python API is quickly evolving.
For regular data analysis, it is recommended to use the scripts which are stable, independently of the underlying Python API.
2D Ptychography (operator-based)#
This is the Ptychography class, which can be used with operators
- class pynx.ptycho.ptycho.OperatorPtycho(lazy=False)#
Base class for an operator on Ptycho2D objects.
- timestamp_increment(p)#
Increment the timestamp counter corresponding to the processing language used (OpenCL or CUDA) Virtual, must be derived.
- Parameters:
w – the object (e.g. wavefront) the operator will be applied to.
- Returns:
- class pynx.ptycho.ptycho.Ptycho(probe=None, obj=None, background=None, data=None, nb_frame_total=None)#
Class for two-dimensional ptychographic data: object, probe, and observed diffraction. This may include only part of the data from a larger dataset.
- Parameters:
probe – the starting estimate of the probe, as a complex 2D numpy array - can be 3D if modes are used. the probe should be centered in the center of the array.
obj – the starting estimate of the object, as a complex 2D numpy array - can be 3D if modes are used.
background – 2D array with the incoherent background.
data – the PtychoData object with all observed frames, ptycho positions
nb_frame_total – total number of frames (used for normalization)
- calc_regularisation_scale()#
Calculate the scale factor for object and probe regularisation. Calculated according to Thibault & Guizar-Sicairos 2012 :return: nothing
- calc_scan_area()#
Compute the scan area for the object and probe, using scipy ConvexHull. The scan area for the object is augmented by twice the average distance between scan positions for a more realistic estimation. scan_area_points is also computed, corresponding to the outline of the scanned area.
- Returns:
Nothing. self.scan_area_probe and self.scan_area_obj are updated, as 2D arrays with the same shape as the object and probe, with False outside the scan area and True inside.
- from_pu()#
Get all relevant arrays from processing unit, if necessary :return: Nothing
- get_background()#
Get the background data array. This will automatically get the latest data, either from GPU or from the host memory, depending where the last changes were made.
- Returns:
the 2D numpy data array
- get_illumination_obj()#
Get the sum of the probe intensity for all illuminations, which is used for object normalisation.
- Returns:
the array of the illumination norm, with the same 2D shape as the object
- get_iobs(i=None, shift=False)#
Get the observed intensity array.
- Parameters:
i – the desired frame. if None, all data is returned.
shift – if True, the array will be shifted so that the array is centered on the frame center (and not the corner).
- Returns:
a masked array of the observed intensity for frame i
- get_llk(noise=None, norm=True)#
Get the log-likelihood.
- Parameters:
noise – noise model, either ‘poisson’, ‘gaussian’ or ‘euclidian’. If None, a dictionary is returned.
norm – if True (the default), the LLK values are normalised
- Returns:
either a single LLK value, or a dictionary
- get_obj(remove_obj_phase_ramp=False)#
Get the object data array. This will automatically get the latest data, either from GPU or from the host memory, depending where the last changes were made.
- Parameters:
remove_obj_phase_ramp – if True, the object will be returned after removing the phase ramp coming from the imperfect centring of the diffraction data (sub-pixel shift). Calculated diffraction patterns using such a corrected object will present a sub-pixel shift relative to the diffraction data. The ramp information comes from the PtychoData phase_ramp_d{x,y} attributes, and should have been calculated beforehand using a ZroPhaseRamp(obj=True) operator.
- Returns:
the 3D numpy data array (nb object modes, nyo, nxo)
- get_obj_coord()#
Get the object coordinates :return: a tuple of two arrays corresponding to the x (columns) and y (rows coordinates)
- get_probe()#
Get the probe data array. This will automatically get the latest data, either from GPU or from the host memory, depending where the last changes were made.
- Returns:
the 3D probe numpy data array
- get_probe_coord()#
Get the probe coordinates :return: a tuple of two arrays corresponding to the x (columns) and y (rows coordinates)
- get_psi(shift=False)#
Get the current Psi array. Note that this is not not intended for an efficient data access, but only for occasional ones.
- Parameters:
shift – if True, the array will be shifted so that the array is centered on the frame center (and not the corner).
- Returns:
the raw psi numpy array for the currently processed stack of frames. Can also be None if no calculation took place.
- get_scan_area_obj()#
Return the scan_area_obj (2D array with the object shape, True inside the area scanned, and False outside). It is computed if necessary. :return: self.scan_area_obj
- get_scan_area_points()#
Return the scan_area_points (outside polygon points delimiting the scan area). It is computed if necessary. :return: self.scan_area_points
- get_scan_area_probe()#
Return the scan_area_probe (2D array with the probe shape, True inside the area scanned, and False outside). It is computed if necessary. :return: self.scan_area_probe
- load_obj_probe_cxi(filename, entry=None, verbose=True)#
Load object and probe from a CXI file, result of a previous optimisation. If no data is already present in the current object, then the pixel size and energy/wavelength are also loaded, and a dummy (one frame) data object is created.
- Parameters:
filename – the CXI filename from which to load the data
entry – the entry to be read. By default, the last in the file is loaded. Can be ‘entry_1’, ‘entry_10’… An integer n can also be supplied, in which case ‘entry_%d’ % n will be read
- Returns:
- print(*args, **kwargs)#
Print function which can be muted e.g. for non-master MPI processes :param args: arguments passed to the print function :param kwargs: keyword arguments passed to the print function :return: nothing
- reset_history()#
Reset history, and set current cycle to zero :return: nothing
- save_obj_probe_cxi(filename, sample_name=None, sample_nx_dict=None, experiment_id=None, instrument=None, note=None, process=None, append=False, shift_phase_zero=False, params=None, remove_obj_phase_ramp=False, cxi_data='object_probe')#
Save the result of the optimisation (object, probe, scan areas) to an HDF5 CXI-like file.
- Parameters:
filename – the file name to save the data to
sample_name – optional, the sample name. Ignored if sample_nx_dict is given.
sample_nx_dict – NXsample object as a dictionary, from silx.io.dictdump.nxtodict.
experiment_id – the string identifying the experiment, e.g.: ‘HC1234: Siemens star calibration tests’
instrument – the string identifying the instrument, e.g.: ‘ESRF id10’
note – a string with a text note giving some additional information about the data, a publication…
process – a dictionary of strings which will be saved in ‘/entry_N/data_1/process_1’. A dictionary entry can also be a ‘note’ as keyword and a dictionary as value - all key/values will then be saved as separate notes. Example: process={‘program’: ‘PyNX’, ‘note’:{‘llk’:1.056, ‘nb_photons’: 1e8}}
append – by default (append=False), any existing file will be overwritten, and the result will be saved as ‘entry_1’. If append==True and the file exists, a new entry_N will be saved instead. This can be used to export the different steps of an optimisation.
shift_phase_zero – if True, remove the linear phase ramp from the object
params – a dictionary of parameters to be saved into process_1/configuration NXcollection
remove_obj_phase_ramp – if True, the object will be saved after removing the phase ramp coming from the imperfect centring of the diffraction data (sub-pixel shift). Calculated diffraction patterns using such a corrected object will present a sub-pixel shift relative to the diffraction data. The ramp information comes from the PtychoData phase_ramp_d{x,y} attributes, and are not re-calculated.
cxi_data –
a string describing what large data arrays are actually saved. This only affects the large object and probe arrays.
’object_probe’: all object and probe modes are saved as complex objects.
’object_phase’: only the first mode of the object’s phase is saved (as float16). This is the format which saves the most space
’object’: save only the complex object
’probe’: save only the complex probe
- Returns:
Nothing. a CXI file is created
- set_background(background)#
Set the incoherent background data array. It will be shifted so that the center of the diffraction image is at (0,0), like the stored intensity.
- Parameters:
background – the background (float32 numpy array)
- Returns:
nothing
- set_obj(obj)#
Set the object data array.
- Parameters:
obj – the object (complex64 numpy array)
- Returns:
nothing
- set_obj_zero_phase_mask(mask)#
Set an object mask, which has the same 2D shape as the object, where values of 1 indicate that the area corresponds to vacuum (or air), and 0 corresponds to some material. Values between 0 and 1 can be given to smooth the transition. This mask will be used to restrain the corresponding area to a null phase, dampening the imaginary part at every object update. :param mask: a floating-point array with the same 2D shape as the object, where values of 1 indicate that the area corresponds to vacuum (or air), and 0 corresponds to the sample. :return: nothing
- set_probe(probe)#
Set the probe data array.
- Parameters:
probe – the probe (complex64 numpy array)
- Returns:
nothing
- update_history(mode='llk', update_obj=False, update_probe=False, update_background=False, update_pos=False, verbose=False, **kwargs)#
Update the history record.
- Parameters:
mode – either ‘llk’ (will record new log-likelihood and number of photons) or ‘algorithm’ (will only update the algorithm) - for the latter case, algorithm should be given as a keyword argument.
verbose – if True, print some info about current process (only if mode==’llk’)
kwargs – other parameters to be recorded, e.g. probe_inertia=xx, algorithm=’DM’
- Returns:
nothing
- class pynx.ptycho.ptycho.PtychoData(iobs=None, positions=None, detector_distance=None, mask=None, pixel_size_detector=None, wavelength=None, near_field=False, padding=0, vidx=None)#
Class for two-dimensional ptychographic data: observed diffraction and probe positions. This may include only part of the data from a larger dataset.
- Parameters:
iobs – 3d array with (nb_frame, ny,nx) shape observed intensity (assumed to follow Poisson statistics). The frames will be stored fft-shifted so that the center of diffraction lies in the (0,0) corner of each image. The supplied frames should have the diffraction center in the middle of the frames. Intensities must be >=0, except for masked pixels which should be <0. For the purpose of Paganin and CTF operators, masked values should be set to inter/extrapolated values and stored as -1-I_interp
positions – (x, y, z) tuple or 2-column array with ptycho probe positions in meters. For 2D data, z is ignored and can be None or missing, e.g. with (x, y) The orthonormal coordinate system follows the CXI/NeXus/McStas convention, with z along the incident beam propagation direction and y vertical towards ceiling. Values for x or y larger than 1e10 indicate that the direct beam is measured, allowing an absolute reference for the probe. Alternatively, those values may also be masked in a numpy masked array.
detector_distance – detector distance in meters
mask – obsolete, should not be used anymore. Masked intensity should be supplied as negative intensities.
pixel_size_detector – in meters, assuming square pixels
wavelength – wavelength of the experiment, in meters.
near_field – True if using near field ptycho
padding – an integer value indicating the number of zero-padded pixels to be used on each side of the observed frames. This can be used for near field ptychography. The input iobs should already padded, the corresponding pixels will be added to the mask.
vidx – array of indices of the positions. This is useful when only a subset of positions is studied, e.g. when the dataset is split between parallel processes.
kwargs – used to detect obsoleted parameters such as mask
- get_required_obj_shape(margin=16)#
Estimate the required object shape
- Parameters:
margin – number of pixels on the border to avoid stepping out of the object, e.g. when refining positions.
- Returns:
(nyo, nxo), the 2D object shape
- pixel_size_object()#
Get the x and y pixel size in object space after a FFT. :return: a tuple (pixel_size_x, pixel_size_y) in meters
- pynx.ptycho.ptycho.algo_string(algo_base, p, update_object, update_probe, update_background=False, update_pos=False)#
Get a short string for the algorithm being run, e.g. ‘DM/o/3p’ for difference map with 1 object and 3 probe modes.
- Parameters:
algo_base – ‘AP’ or ‘ML’ or ‘DM’
p – the ptycho object
update_object – True if updating the object
update_probe – True if updating the probe
update_background – True if updating the background
update_pos – True if updating the positions
- Returns:
a short string for the algorithm
- pynx.ptycho.ptycho.calc_throughput(p: Ptycho = None, cxi=None, verbose=False)#
Analyse the throughput after a series of algorithms, either from a Ptycho object or from a CXI file. A few things like object & probe smoothing are not taken into account, or the use of bilinear interpolation, background… :param p: the Ptycho object the timings are extracted from. :param cxi: the CXI file the history of cycles will be obtained from :param verbose: if True, print the average throughput per algorithm step :return: the average throughput in Gbyte/s
- pynx.ptycho.ptycho.save_obj_probe_cxi(filename, obj, probe, wavelength, detector_distance, pixel_size_detector, llk_poisson, llk_gaussian, llk_euclidian, nb_photons_calc, history, pixel_size_object, posxy, posxy_c, posxy0, scale=None, obj_zero_phase_mask=None, scan_area_obj=None, scan_area_probe=None, obj_illumination=None, background=None, sample_name=None, sample_nx_dict=None, experiment_id=None, instrument=None, note=None, process=None, append=False, shift_phase_zero=False, params=None, extra_data=None, obj_phase_ramp=None)#
Save the result of the optimisation (object, probe, scan areas) to an HDF5 CXI-like file.
Note that object and probed are flipped (up/down) to have a (top, left) array origin.
- Parameters:
filename – the file name to save the data to
obj – the object to save
probe – the probe to save
wavelength – the wavelength (SI unit)
detector_distance – detector distance
pixel_size_detector – the detector’s pixel size
llk_euclidian (llk_poisson, llk_gaussian,) – normalised log-likelihood values
pixel_size_object – the object pixel size
posxy – the final scanning positions (SI units)
posxy0 – the initial scanning positions (SI units)
posxy_c – xy coordinates of the object center
scale – array of per-frame scaling factir
obj_zero_phase_mask – the area used for the object phase optimisation (maybe obsolete..)
scan_area_obj – scan area mask on the object array
scan_area_probe – scan area mask on the probe array
obj_illumination – integrated incident intensity on the object area
background – incoherent background
sample_name – optional, the sample name. Ignored if sample_nx_dict is given.
sample_nx_dict – NXsample object as a dictionary, from silx.io.dictdump.nxtodict.
experiment_id – the string identifying the experiment, e.g.: ‘HC1234: Siemens star calibration tests’
instrument – the string identifying the instrument, e.g.: ‘ESRF id10’
note – a string with a text note giving some additional information about the data, a publication…
process – a dictionary of strings which will be saved in ‘/entry_N/data_1/process_1’. A dictionary entry can also be a ‘note’ as keyword and a dictionary as value - all key/values will then be saved as separate notes. Example: process={‘program’: ‘PyNX’, ‘note’:{‘llk’:1.056, ‘nb_photons’: 1e8}}
append – by default (append=False), any existing file will be overwritten, and the result will be saved as ‘entry_1’. If append==True and the file exists, a new entry_N will be saved instead. This can be used to export the different steps of an optimisation.
shift_phase_zero – if True, centre the object phase around zero
params – a dictionary of parameters to be saved into process_1/configuration NXcollection
extra_data – a dictionary of data which will be saved as entry/_extra_data, and may be useful for debugging. each value may itself be a dictionary of values to save
obj_phase_ramp – the shifts (dx, dy) of the average calculated intensity from the array centre
- Returns:
Nothing. a CXI file is created
- pynx.ptycho.ptycho.save_obj_probe_nxtomo(filename, obj, probe, wavelength, detector_distance, pixel_size_detector, llk_poisson, llk_gaussian, llk_euclidian, nb_photons_calc, history, tomo_angle, crop_shape, pixel_size_object, posxy, posxy_c, posxy0, scale=None, obj_zero_phase_mask=None, scan_area_obj=None, scan_area_probe=None, obj_illumination=None, background=None, sample_name=None, experiment_id=None, instrument=None, note=None, process=None, append=False, shift_phase_zero=False, params=None, extra_data=None, obj_phase_ramp=None, deriv_phase=True, retry=True)#
Save the result of the optimisation (object, probe, scan areas) to an NXTomo (https://manual.nexusformat.org/classes/applications/NXtomo.html) file. This is used for ptycho-tomo experiments, with each projection (tomo angle) being written independently and possibly concurrently.
Note that object and probed are flipped (up/down) to have a (top, left) array origin.
- Parameters:
filename – the file name to save the data to
obj – the object to save
probe – the probe to save (will not be saved if None)
wavelength – the wavelength (SI unit)
detector_distance – detector distance
pixel_size_detector – the detector’s pixel size
llk_euclidian (llk_poisson, llk_gaussian,) – normalised log-likelihood values
nb_photons_calc – total number of photons, calculated
history – the pynx.utils.history.History record of the optimisation
tomo_angle – the rotation angle for this projection, in radians
crop_shape – the 2D shape (ny, nx) the object should be cropped to before saving, so all projections have the same shape.
pixel_size_object – the object pixel size
posxy – the final scanning positions (SI units)
posxy0 – the initial scanning positions (SI units)
posxy_c – xy coordinates of the object center
background – incoherent background
sample_name – the sample name. If not given, the filename without extension will be used.
experiment_id – the string identifying the experiment, e.g.: ‘HC1234: Siemens star calibration tests’
instrument – the string identifying the instrument, e.g.: ‘ESRF id10’
note – a string with a text note giving some additional information about the data, a publication…
process – a dictionary of strings which will be saved in ‘/entry_N/data_1/process_1’. A dictionary entry can also be a ‘note’ as keyword and a dictionary as value - all key/values will then be saved as separate notes. Example: process={‘program’: ‘PyNX’, ‘note’:{‘llk’:1.056, ‘nb_photons’: 1e8}}
params – a dictionary of parameters to be saved into process_1/configuration NXcollection
extra_data – a dictionary of data which will be saved as entry/_extra_data, and may be useful for debugging. each value may itself be a dictionary of values to save
obj_phase_ramp – the shifts (dx, dy) of the average calculated intensity from the array centre
deriv_phase – if True, the phase derivative (along the horizontal direction) will be computed from the complex object.
retry – if True, the opening of the .nx file in append mode will be re-tried until successful - this is needed for independent, asynchronous parallel optimisations.
- Returns:
Nothing. A nexus file is created or appended
- pynx.ptycho.ptycho.save_ptycho_data_cxi(file_name, iobs, pixel_size, wavelength, detector_distance, x, y, z=None, monitor=None, mask=None, dark=None, instrument='', overwrite=False, scan=None, params=None, verbose=False, **kwargs)#
Save the Ptychography scan data using the CXI format (see http://cxidb.org)
- Parameters:
file_name – the file name (including relative or full path) to save the data to
iobs – the observed intensity, with shape (nb_frame, ny, nx)
pixel_size – the detector pixel size in meters
wavelength – the experiment wavelength
x – the x scan positions
y – the y scan positions
z – the z scan positions (default=None)
monitor – the monitor
mask – the mask for the observed frames
dark – the incoherent background
instrument – a string with the name of the instrument (e.g. ‘ESRF id16A’)
overwrite – if True, will overwrite an existing file
params – a dictionary of parameters which will be saved as a NXcollection
verbose – if True, print something.
- Returns:
2D Ptychography Operators#
This section lists the operators, which can be imported automatically using from pynx.ptycho import * or from pynx.ptycho.operator import *. When this import is done, either CUDA (preferably) or OpenCL operators will be imported. The documentation below corresponds to OpenCL operators, but this should be identical to CUDA operators.
- class pynx.ptycho.cl_operator.AP(update_object=True, update_probe=False, update_background=False, floating_intensity=False, nb_cycle=1, calc_llk=False, show_obj_probe=False, fig_num=-1, obj_smooth_sigma=0, obj_inertia=0.01, probe_smooth_sigma=0, probe_inertia=0.001, update_pos=False, pos_mult=1, pos_max_shift=2, pos_min_shift=0, pos_threshold=0.2, pos_history=False, zero_phase_ramp=True, background_smooth_sigma=0)#
Perform a complete Alternating Projection cycle: - forward all object*probe views to Fourier space and apply the observed amplitude - back-project to object space and project onto (probe, object) - update background optionally
- Parameters:
update_object – update object ?
update_probe – update probe ?
update_background – update background ?
floating_intensity – optimise floating intensity scale factor
nb_cycle – number of cycles to perform. Equivalent to AP(…)**nb_cycle
calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle
show_obj_probe – if a positive integer number N, the object & probe will be displayed every N cycle. By default 0 (no plot)
fig_num – the number of the figure to plot the object and probe, as for ShowObjProbe()
obj_smooth_sigma – if > 0, the previous object array (used for inertia) will convoluted by a gaussian array of this standard deviation.
obj_inertia – the updated object retains this relative amount of the previous object.
probe_smooth_sigma – if > 0, the previous probe array (used for inertia) will convoluted by a gaussian array of this standard deviation.
probe_inertia – the updated probe retains this relative amount of the previous probe.
update_pos – if True, update positions. This automatically inhibits object and probe update
pos_max_shift – maximum allowed shift (in pixels) per scan position (default=2)
pos_min_shift – minimum required shift (in pixels) per scan position (default=0)
pos_threshold – if the integrated grad_obj*probe along dx or dy is smaller than (grad_obj*probe).mean()*threshold, then the shift is ignored. This allows to prevent position shifts where there is little contrast in the object. It applies independently to x and y.
pos_mult – multiply the calculated position shifts by this value. Useful since the calculated shifts usually are a fraction of the actual shift.
pos_history – if True, save the position history (for debugging)
zero_phase_ramp – if True, the conjugate phase ramp in the object and probe will be removed by centring the FT of the probe, at the end and before every display. Ignored for near field.
background_smooth_sigma – gaussian convolution parameter for the background update
- class pynx.ptycho.cl_operator.ApplyAmplitude(calc_llk=False, update_background=False, scale_in=1, scale_out=1, sum_icalc=False, background_smooth_sigma=0)#
Apply the magnitude from observed intensities, keep the phase. Applies to a stack of N views.
- Parameters:
calc_llk – if True, the log-likelihood will be calculated for this stack.
update_background – if True, update the background. The new background is automatically updated once the last stack is processed.
scale_in – a scale factor by which the input values should be multiplied, typically because of FFT
scale_out – a scale factor by which the output values should be multiplied, typically because of FFT
sum_icalc – if True, will sum the calculated intensity on each frame, and store this in p._cl_icalc_sum, to be used for floating intensities
background_smooth_sigma – sigma for gaussian smoothing of the background
- class pynx.ptycho.cl_operator.ApplyPhaseRamp(dx, dy, probe=False, obj=False)#
Apply a given phase ramp to the object and/or probe. The actual phase factor is: For the probe: np.exp(-2j*np.pi*(x * dx + y*dy)) For the object: np.exp(2j*np.pi*(x * dx + y*dy)) Where (x,y) are reduced pixel coordinates going from -.5 to .5 for the probe, and +/-0.5 * nxo/nx (or nyo/ny) for the object.
If the phase ramp is applied to both object and probe, the calculated intensity patterns remain unchanged.
- Parameters:
dy (dx,) – relative shifts from the centre, calculated in reciprocal space (probe array pixel coordinates)
obj – if True, apply the correction to the object
probe – if True, apply the correction to the probe
- class pynx.ptycho.cl_operator.BackgroundFilter(sigma)#
Apply a Gaussian filter to the background array. This operator will perform a FT on the backround, multiply it by a Gaussian, and back-FT the array. The resulting array is normalised so that its sum is equal to the original array sum.
- Parameters:
sigma – standard deviation for the Gaussian kernel in the detector space, corresponding to a convolution of the background by a Gaussian of FWHM=2.3548*sigma
- class pynx.ptycho.cl_operator.CLObsDataStack(cl_obs, cl_x, cl_y, i, npsi)#
Class to store a stack (e.g. 16 frames) of observed data in OpenCL space
- Parameters:
cl_obs – pyopencl array of observed data, with N frames
cl_y (cl_x,) – pyopencl arrays of the positions (in pixels) of the different frames
i – index of the first frame
npsi – number of valid frames (others are filled with zeros)
- class pynx.ptycho.cl_operator.CLOperatorPtycho(processing_unit=None)#
Base class for a operators on Ptycho objects using OpenCL
- apply_ops_mul(pty)#
Apply the series of operators stored in self.ops to a Ptycho2D object. The operators are applied one after the other.
- Parameters:
pty – the Ptycho2D to which the operators will be applied.
- Returns:
the Ptycho2D object, after application of all the operators in sequence
- init_cl_vobs(p)#
Initialize observed intensity and scan positions in OpenCL space
- Parameters:
p – the Ptycho object the operator will be applied to.
- Returns:
- prepare_data(p: Ptycho)#
Make sure the data to be used is in the correct memory (host or GPU) for the operator. Virtual, must be derived.
- Parameters:
p – the Ptycho object the operator will be applied to.
- Returns:
- timestamp_increment(p)#
Increment the timestamp counter corresponding to the processing language used (OpenCL or CUDA) Virtual, must be derived.
- Parameters:
w – the object (e.g. wavefront) the operator will be applied to.
- Returns:
- view_copy(pty, i_source, i_dest)#
Create a new view of the object by copying the original data. This will make a copy of all relevant data, which can be a wavefront, CDI object, Ptychography object, probe and psi arrais, etc…
This (virtual) function is used to make temporary copies of the data Operators apply to. This is used to make linear combination (additions) of operators, which requires several copies of the data. As the copying part depends on the processing unit used (GPU, CPU) and the exact data to duplicate, this is a pure virtual class, which must be derived to be used. .. note:
- this should only copy the 'active' data (which is affected by calculations) - index 0 corresponds to the original array, to which subsequent operators will be applied to
- Parameters:
obj – the object where the data will be duplicated
i_source – the index (integer) of the source object data
i_dest – the index (integer) of the destination object data
- Returns:
nothing. The object is modified in-place
- view_purge(pty, i)#
Purge the different views of an object (except the main one). As the purging depends on the processing unit used (GPU, CPU) and the exact data to purge, this is a pure virtual function, which must be derived to be used. :param obj: the object where the view will be purged :param i_view: the index of the view to purge. If None, all views are purged. This de-registers the view. :return: nothing
- view_register(obj)#
Creates a new unique view key in an object. When finished with this view, it should be de-registered using view_purge. Note that it only reserves the key, but does not create the view. :return: an integer value, which corresponds to yet-unused key in the object’s view.
- view_sum(pty, i_source, i_dest)#
Add the view data from one index into another. As the summing depends on the processing unit used (GPU, CPU) and the exact data to sum, this is a pure virtual function, which must be derived to be used. :param obj: the object where the data will be duplicated :param i_source: the index (integer) of the source object data :param i_dest: the index (integer) of the destination object data :return: nothing. The object is modified in-place
- view_swap(pty, i1, i2)#
Swap the object view between index i1 and i2. As the swapping part depends on the processing unit used (GPU, CPU) and the exact data to swap, this is a pure virtual function, which must be derived to be used. :param obj: the object where the data will be duplicated :param i1: the index (integer) of the first object data :param i2: the index (integer) of the second object data :return: nothing. The object is modified in-place
- class pynx.ptycho.cl_operator.CLOperatorPtychoPower(op, n)#
- class pynx.ptycho.cl_operator.CLOperatorPtychoSum(op1, op2)#
- class pynx.ptycho.cl_operator.CLProcessingUnitPtycho#
Processing unit in OpenCL space, for 2D Ptycho operations.
Handles initializing the context and kernels.
- cl_init_kernel_n(n)#
Init kernels specifically written for reduction with N-sized array (N being normally the stack_size) :param n: the size for the float_n arrays used :return: Nothing.
- cl_init_kernels()#
Initialize opencl kernels :return: nothing
- set_stack_size(s)#
Change the number of frames which are stacked to perform all operations in //. If it is larger than the total number of frames, operators like AP, DM, ML will loop over all the stacks. :param s: an integer number (default=16) :return: nothing
- class pynx.ptycho.cl_operator.Calc2Obs(nb_photons_per_frame=None, poisson_noise=False, scale_probe_direct_beam=True)#
Copy the calculated intensities to the observed ones, optionally with Poisson noise. This operator will loop other all stacks of frames, multiply object and probe and propagate the wavefront to the detector, and compute the new intensities.
- Parameters:
nb_photons_per_frame – average number of photons per frame, to scale the images. If None, no scaling is performed.
poisson_noise – if True, Poisson noise will be applied on the calculated frames. This uses numpy.random.poisson and is not GPU-optimised.
scale_probe_direct_beam – if True, poisson_noise is True, and the dataset includes at least one direct beam frame (without sample), then the probe will also be scaled so that it corresponds to the direct beam intensity.
- op(p: Ptycho)#
Applies the operator to one object. Virtual function, must be derived. By default this is the identity operator.
- Returns:
the result of the operation, usually the same type as the input, which is modified in-place. But this can also be a scalar (reduction).
- timestamp_increment(p)#
Increment the timestamp counter corresponding to the processing language used (OpenCL or CUDA) Virtual, must be derived.
- Parameters:
w – the object (e.g. wavefront) the operator will be applied to.
- Returns:
- class pynx.ptycho.cl_operator.Calc2Obs1#
Copy the calculated intensities to the observed ones. Can be used for simulation. Applies to a stack of N views, assumes the current Psi is already in Fourier space.
- class pynx.ptycho.cl_operator.CalcIllumination(processing_unit=None)#
Compute the integrated illumination of the object by all probe positions
- op(p: Ptycho)#
- Parameters:
p – the Ptycho object this operator applies to
- Returns:
the updated Ptycho object
- timestamp_increment(p)#
Increment the timestamp counter corresponding to the processing language used (OpenCL or CUDA) Virtual, must be derived.
- Parameters:
w – the object (e.g. wavefront) the operator will be applied to.
- Returns:
- class pynx.ptycho.cl_operator.CenterObjProbe(max_shift=5, power=2, verbose=False)#
Operator to check the center of mass of the probe and shift both object and probe if necessary.
- Parameters:
max_shift – the maximum shift of the probe with respect to the center of the array, in pixels. The probe and object are only translated if the shift is larger than this value.
power – the center of mass is calculated on the amplitude of the array elevated at this power.
verbose – print deviation if verbose=True
- class pynx.ptycho.cl_operator.DM(update_object=True, update_probe=True, update_background=False, nb_cycle=1, calc_llk=False, show_obj_probe=False, fig_num=-1, obj_smooth_sigma=0, obj_inertia=0.01, probe_smooth_sigma=0, probe_inertia=0.001, center_probe_n=0, center_probe_max_shift=5, loop_obj_probe=1, update_pos=False, pos_mult=1, pos_max_shift=2, pos_min_shift=0, pos_threshold=0.2, pos_history=False, zero_phase_ramp=True, background_smooth_sigma=0, alpha=0.02)#
Operator to perform a complete Difference Map cycle, updating the Psi views for all stack of frames, as well as updating the object and/or probe.
- Parameters:
update_object – update object ?
update_probe – update probe ?
update_background – update background ?
nb_cycle – number of cycles to perform. Equivalent to DM(…)**nb_cycle
calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle
show_obj_probe – if a positive integer number N, the object & probe will be displayed every N cycle. By default 0 (no plot)
fig_num – the number of the figure to plot the object and probe, as for ShowObjProbe()
obj_smooth_sigma – if > 0, the previous object array (used for inertia) will convoluted by a gaussian array of this standard deviation.
obj_inertia – the updated object retains this relative amount of the previous object.
probe_smooth_sigma – if > 0, the previous probe array (used for inertia) will convoluted by a gaussian array of this standard deviation.
probe_inertia – the updated probe retains this relative amount of the previous probe.
center_probe_n – test the probe every N cycle for deviation from the center. If deviation is larger than center_probe_max_shift, probe and object are shifted to correct. If 0 (the default), the probe centering is never calculated.
center_probe_max_shift – maximum deviation from the center (in pixels) to trigger a position correction
loop_obj_probe – when both object and probe are updated, it can be more stable to loop the object and probe update for a more stable optimisation, but slower.
update_pos – positive integer, if >0, update positions every ‘update_pos’ cycle. Note that object and probe are not updated during the same cycle as positions. Still experimental, recommended are 5, 10 (default=False, positions are not updated).
pos_max_shift – maximum allowed shift (in pixels) per scan position (default=2)
pos_min_shift – minimum required shift (in pixels) per scan position (default=0)
pos_threshold – if the integrated grad_obj*probe along dx or dy is smaller than (grad_obj*probe).mean()*threshold, then the shift is ignored. This allows to prevent position shifts where there is little contrast in the object. It applies independently to x and y.
pos_mult – multiply the calculated position shifts by this value. Useful since the calculated shifts usually are a fraction of the actual shift.
pos_history – if True, save the position history (for debugging)
zero_phase_ramp – if True, the conjugate phase ramp in the object and probe will be removed by centring the FT of the probe, at the end and before every display. Ignored for near field.
background_smooth_sigma – gaussian convolution parameter for the background update
alpha – mixing parameter between AP and DM, to bring more stability to DM. Default is 0.02, equivalent to 0.02*AP+0.98*DM.
- class pynx.ptycho.cl_operator.DM1(alpha=0.0)#
Equivalent to operator: (2-alpha) * ObjProbe2Psi() - (1-alpha)
- Parameters:
alpha – alpha parameter to mix some AP within DM
- class pynx.ptycho.cl_operator.DM2(alpha=0.0)#
# Psi(n+1) = (1-alpha) * (Psi(n) - P*O) + Psi_fourier
This operator assumes that Psi_fourier is the current Psi, and that Psi(n) is in p._cl_psi_copy
On output Psi(n+1) is the current Psi
- Parameters:
alpha – alpha parameter to mix AP with DM
- class pynx.ptycho.cl_operator.FT(scale=True, copy_psi=False)#
Forward Fourier-transform a Psi array, i.e. a stack of N Obj*Probe views
- Parameters:
scale – if True, the FFT will be normalized.
copy_psi – if True, perform the FT out-of-place, keeping the original Psi array in the Ptycho object’s _cl_psi_copy
- class pynx.ptycho.cl_operator.FreePU(processing_unit=None)#
Operator freeing OpenCL memory. The fft plan/app in self.processing_unit is removed, as well as any OpenCL pyopencl.array.Array attribute in the supplied wavefront.
- op(p: Ptycho)#
- Parameters:
p – the Ptycho object this operator applies to
- Returns:
the updated Ptycho object
- timestamp_increment(p)#
Increment the timestamp counter corresponding to the processing language used (OpenCL or CUDA) Virtual, must be derived.
- Parameters:
w – the object (e.g. wavefront) the operator will be applied to.
- Returns:
- class pynx.ptycho.cl_operator.Grad(update_object=True, update_probe=False, update_background=False, floating_intensity=False, reg_fac_obj=0, reg_fac_probe=0, calc_llk=False)#
Operator to compute the object and/or probe and/or background gradient. The gradient is stored in the ptycho object. It is assumed that the GPU gradient arrays have been already created, normally by the calling ML operator.
- Parameters:
update_object – compute gradient for the object ?
update_probe – compute gradient for the probe ?
update_background – compute gradient for the background ?
floating_intensity – optimise floating intensity scale factor
calc_llk – calculate llk while in Fourier space
- class pynx.ptycho.cl_operator.IFT(scale=True)#
Backward Fourier-transform a Psi array, i.e. a stack of N Obj*Probe views
- Parameters:
scale – if True, the FFT will be normalized.
- class pynx.ptycho.cl_operator.LLK(scale=False)#
Log-likelihood reduction kernel. Can only be used when Psi is in diffraction space. This is a reduction operator - it will write llk as an argument in the Ptycho object, and return the object. If _cl_stack_i==0, the llk is re-initialized. Otherwise it is added to the current value.
The LLK can be calculated directly from object and probe using: p = LoopStack(LLK() * FT() * ObjProbe2Psi()) * p
- Parameters:
scale – if True, will scale the calculated amplitude to calculate the log-likelihood. The amplitudes are left unchanged.
- class pynx.ptycho.cl_operator.LoopScan(op)#
Operator to apply a given operator sequentially to different scans,
- Parameters:
op – the operator to apply, which can be a multiplication of operators
- class pynx.ptycho.cl_operator.LoopStack(op, keep_psi=False, keep_psi_copy=False)#
Operator to apply a given operator sequentially to the complete stack of frames of a ptycho object.
Make sure that the current selected stack is in a correct state (i.e. in sample or detector space,…) before starting such a loop with keep_psi=True.
- Parameters:
op – the operator to apply, which can be a multiplication of operators
keep_psi – if True, when switching between stacks, store psi in p._cl_psi_v.
keep_psi_copy – if True, when switching between stacks, store and restore psi_copy in p._cl_psi_copy_v. Used for algorithms (RAAR) requiring two psi copies across cycles.
- class pynx.ptycho.cl_operator.ML(nb_cycle=1, update_object=True, update_probe=False, update_background=False, floating_intensity=False, reg_fac_obj=0, reg_fac_probe=0, calc_llk=False, show_obj_probe=False, fig_num=-1, update_pos=False, pos_mult=1, pos_max_shift=2, pos_min_shift=0, pos_threshold=0.2, pos_history=False, zero_phase_ramp=True)#
Operator to perform a maximum-likelihood conjugate-gradient minimization.
- Parameters:
update_object – update object ?
update_probe – update probe ?
update_background – update background ?
floating_intensity – optimise floating intensity scale factor
reg_fac_obj – use this regularization factor for the object (if 0, no regularization)
reg_fac_probe – use this regularization factor for the probe (if 0, no regularization)
calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle
show_obj_probe – if a positive integer number N, the object & probe will be displayed every N cycle. By default 0 (no plot)
fig_num – the number of the figure to plot the object and probe, as for ShowObjProbe()
update_pos – positive integer, if >0, update positions every ‘update_pos’ cycle. Note that object and probe are not updated during the same cycle as positions. Still experimental, recommended are 5, 10 (default=False, positions are not updated).
pos_max_shift – maximum allowed shift (in pixels) per scan position (default=2)
pos_min_shift – minimum required shift (in pixels) per scan position (default=0)
pos_threshold – if the integrated grad_obj*probe along dx or dy is smaller than (grad_obj*probe).mean()*threshold, then the shift is ignored. This allows to prevent position shifts where there is little contrast in the object. It applies independently to x and y.
pos_mult – multiply the calculated position shifts by this value. Useful since the calculated shifts usually are a fraction of the actual shift.
pos_history – if True, save the position history (for debugging)
zero_phase_ramp – if True, the conjugate phase ramp in the object and probe will be removed by centring the FT of the probe, at the beginning and end. Ignored for near field.
- class pynx.ptycho.cl_operator.MemUsage(verbose=True)#
Print memory usage of current process (RSS on host) and used GPU memory
- op(p: Ptycho)#
- Parameters:
p – the ptycho object this operator applies to
- Returns:
the updated ptycho object
- prepare_data(p: Ptycho)#
Make sure the data to be used is in the correct memory (host or GPU) for the operator. Virtual, must be derived.
- Parameters:
p – the Ptycho object the operator will be applied to.
- Returns:
- timestamp_increment(p)#
Increment the timestamp counter corresponding to the processing language used (OpenCL or CUDA) Virtual, must be derived.
- Parameters:
w – the object (e.g. wavefront) the operator will be applied to.
- Returns:
- class pynx.ptycho.cl_operator.ObjProbe2LLK(use_copy=False)#
Compute the LLK from the current object and probe. This includes the object and probe multiplication and propagation
- Parameters:
use_copy – if True, will work on the ptycho object’s _cl_psi_copy and leave _cl_psi untouched
- class pynx.ptycho.cl_operator.ObjProbe2Psi(processing_unit=None)#
Computes Psi = Obj(r) * Probe(r-r_j) for a stack of N probe positions.
- class pynx.ptycho.cl_operator.PropagateApplyAmplitude(calc_llk=False, update_background=False, background_smooth_sigma=0, sum_icalc=False, copy_psi=False)#
Propagate to the detector plane (either in far or near field, perform the magnitude projection, and propagate back to the object plane. This applies to a stack of frames.
- Parameters:
calc_llk – if True, calculate llk while in the detector plane.
update_background – if >0, update the background with the difference between the observed and calculated intensity (with a damping factor), averaged over all frames.
background_smooth_sigma – sigma for the gaussian smoothing of the updated background
sum_icalc – if True, update the background.
copy_psi – if True, perform the propagation out-of-place, keeping the original Psi array in the Ptycho object’s _cl_psi_copy
- class pynx.ptycho.cl_operator.PropagateNearField(forward=True, copy_psi=False)#
Near field propagator
- Parameters:
forward – if True, propagate forward, otherwise backward. The distance is taken from the ptycho data this operator applies to.
copy_psi – if True, perform the propagation out-of-place, keeping the original Psi array in the Ptycho object’s _cl_psi_copy
- class pynx.ptycho.cl_operator.Psi2Obj(floating_intensity=False)#
Computes updated Obj(r) contributions from Psi and Probe(r-r_j), for a stack of N probe positions.
- Parameters:
floating_intensity – if True, the intensity will be considered ‘floating’ from frame to frame, and a scale factor will be adjusted for each frame.
- class pynx.ptycho.cl_operator.Psi2ObjMerge(inertia=0.01, smooth_sigma=0, floating_intensity=False)#
Call this when all stack of probe positions have been processed, and the final update of the object can be calculated. Temporary arrays are cleaned up
- Parameters:
inertia – a regularisation factor to set the object inertia.
floating_intensity – if True, the floating scale factors will be corrected so that the average is 1.
smooth_sigma – if > 0, the previous object array (used for inertia) will be convolved by a gaussian with this sigma.
- class pynx.ptycho.cl_operator.Psi2PosMerge(multiplier=1, max_displ=2, min_displ=0, threshold=0.0, save_position_history=False)#
Merge scan position shifts, once the entire stqck of frames has been processed.
- Parameters:
multiplier – the computed displacements are multiplied by this value
max_displ – the displacements (at each iteration) are capped to this value (in pixels), after applying the multiplier.
min_displ – the displacements (at each iteration) are ignored if smaller than this value (in pixels), after applying the multiplier.
threshold – if the integrated grad_obj*probe along dx or dy is smaller than (grad_obj*probe).mean()*threshold, then the shift is ignored. This allows to prevent shifts where there is little contrast in the object. It applies independently to x and y.
save_position_history – if True, save the position history in the ptycho object (slow, for debugging)
- class pynx.ptycho.cl_operator.Psi2PosShift(processing_unit=None)#
Computes scan position shifts, by comparing the updated Psi array to object*probe, for a stack of frames.
- class pynx.ptycho.cl_operator.Psi2Probe(processing_unit=None)#
Computes updated Probe contributions from Psi and Obj, for a stack of N probe positions.
- class pynx.ptycho.cl_operator.Psi2ProbeMerge(inertia=0.001, smooth_sigma=0)#
Call this when all stack of probe positions have been processed, and the final update of the probe can be calculated. Temporary arrays are cleaned up.
- Parameters:
inertia – a regularisation factor to set the object inertia.
smooth_sigma – if > 0, the previous object array (used for inertia) will be convolved by a gaussian with this sigma.
- class pynx.ptycho.cl_operator.PurgeStacks(processing_unit=None)#
Operator to delete stored psi stacks in a Ptycho object’s _cl_psi_v.
This should be called for each main operator using LoopStack(), once it is finished processing, in order to avoid having another operator using the stored stacks, and to free memory.
- class pynx.ptycho.cl_operator.QuadraticPhase(factor, scale=1)#
Operator applying a quadratic phase factor
Application of a quadratic phase factor, and optionally a scale factor.
The actual factor is: \(scale * e^{i * factor * ((ix/nx)^2 + (iy/ny)^2)}\) where ix and iy are the integer indices of the pixels
- Parameters:
factor – the factor for the phase calculation.
scale – the data will be scaled by this factor. Useful to normalize before/after a Fourier transform, without accessing twice the array data.
- class pynx.ptycho.cl_operator.RAAR(update_object=True, update_probe=True, update_background=False, beta=0.9, nb_cycle=1, calc_llk=False, show_obj_probe=False, fig_num=-1, obj_smooth_sigma=0, obj_inertia=0.01, probe_smooth_sigma=0, probe_inertia=0.001, center_probe_n=0, center_probe_max_shift=5, loop_obj_probe=1, update_pos=False, pos_mult=1, pos_max_shift=2, pos_min_shift=0, pos_threshold=0.2, pos_history=False, zero_phase_ramp=True, background_smooth_sigma=0)#
Operator to perform a complete Relaxed Averaged Alternating Reflectors (RAAR) cycle, updating the Psi views for all stack of frames, as well as the object and/or probe.
- Parameters:
update_object – update object ?
update_probe – update probe ?
update_background – update background ?
beta – RAAR beta parameter.
nb_cycle – number of cycles to perform. Equivalent to DM(…)**nb_cycle
calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle
show_obj_probe – if a positive integer number N, the object & probe will be displayed every N cycle. By default 0 (no plot)
fig_num – the number of the figure to plot the object and probe, as for ShowObjProbe()
obj_smooth_sigma – if > 0, the previous object array (used for inertia) will convoluted by a gaussian array of this standard deviation.
obj_inertia – the updated object retains this relative amount of the previous object.
probe_smooth_sigma – if > 0, the previous probe array (used for inertia) will convoluted by a gaussian array of this standard deviation.
probe_inertia – the updated probe retains this relative amount of the previous probe.
center_probe_n – test the probe every N cycle for deviation from the center. If deviation is larger than center_probe_max_shift, probe and object are shifted to correct. If 0 (the default), the probe centering is never calculated.
center_probe_max_shift – maximum deviation from the center (in pixels) to trigger a position correction
loop_obj_probe – when both object and probe are updated, it can be more stable to loop the object and probe update for a more stable optimisation, but slower.
update_pos – positive integer, if >0, update positions every ‘update_pos’ cycle. Note that object and probe are not updated during the same cycle as positions. Still experimental, recommended are 5, 10 (default=False, positions are not updated).
pos_max_shift – maximum allowed shift (in pixels) per scan position (default=2)
pos_min_shift – minimum required shift (in pixels) per scan position (default=0)
pos_threshold – if the integrated grad_obj*probe along dx or dy is smaller than (grad_obj*probe).mean()*threshold, then the shift is ignored. This allows to prevent position shifts where there is little contrast in the object. It applies independently to x and y.
pos_mult – multiply the calculated position shifts by this value. Useful since the calculated shifts usually are a fraction of the actual shift.
pos_history – if True, save the position history (for debugging)
zero_phase_ramp – if True, the conjugate phase ramp in the object and probe will be removed by centring the FT of the probe, at the end and before every display. Ignored for near field.
background_smooth_sigma – gaussian convolution parameter for the background update
- class pynx.ptycho.cl_operator.RAARProj(beta=0.9)#
# Psi(n+1) = beta/2 [2PO-Psi(n)+ Psi_Fourier] +(1-beta)Psi_Fourier
This operator assumes that Psi_fourier is the current Psi, and that Psi(n) is in p._cl_psi_copy
On output Psi(n+1) is the current Psi
- Parameters:
alpha – beta RAAR parameter
- class pynx.ptycho.cl_operator.Scale(x, obj=True, probe=True, psi=True)#
Multiply the ptycho object by a scalar (real or complex).
- Parameters:
x – the scaling factor
obj – if True, scale the object
probe – if True, scale the probe
psi – if True, scale the all the psi arrays, _cl_psi as well as _cl_psi_v
- class pynx.ptycho.cl_operator.ScaleObjProbe(verbose=False, absolute=True)#
Operator to scale the object and probe so that they have the same magnitude, and that the product of object*probe matches the observed intensity (i.e. sum(abs(obj*probe)**2) = sum(iobs))
- Parameters:
verbose – print deviation if verbose=True
absolute – if True, the absolute scale is computed by comparing the calculated intensities with the observed ones. If False, only the relative scale of object and probe is modified (to avoid numerical divergence), not affecting the calculated intensities.
- class pynx.ptycho.cl_operator.SelectStack(stack_i, keep_psi=False, keep_psi_copy=False)#
Operator to select a stack of observed frames to work on. Note that once this operation has been applied, the new Psi value may be undefined (empty array), if no previous array existed.
Select a new stack of frames, swapping data to store the last calculated psi array in the corresponding, ptycho object’s _cl_psi_v[i] dictionary.
What happens is: * keep_psi=False: only the stack index in p is changed (p._cl_stack_i=stack_i)
- keep_psi=True: the previous psi is stored in p._cl_psi_v[p._cl_stack_i], the new psi is swapped
with p._cl_psi_v[stack_i] if it exists, otherwise initialized as an empty array.
- Parameters:
stack_i – the stack index.
keep_psi – if True, when switching between stacks, store and restore psi in p._cl_psi_v.
keep_psi_copy – if True, when switching between stacks, store and restore psi_copy in p._cl_psi_copy_v. Used for algorithms (RAAR) requiring two psi copies across cycles.
- class pynx.ptycho.cl_operator.SumIntensity1(icalc=None, iobs=None)#
Operator to compute the sum of calculated and/or observed frames. When calculating the sum of observed intensities, masked values will be replaced by calculated ones.
This operator applies to a single stack, and will perform the object*probe multiplication, propagate to the detector, and compute the sums.
- Parameters:
icalc – the GPU array holding the sum of calculated intensities. If None, it is not calculated.
iobs – the GPU array holding the sum of observed intensities. If None, it is not calculated.
- class pynx.ptycho.cl_operator.ZeroPhaseRamp(obj=False)#
Operator to remove the linear phase ramp on both the object and the probe, for far field ptycho. This first computes the center of mass of the square norm of Fourier transform of the probe, and then corrects both probe and object for the phase ramp corresponding to the shift relative to the center of the array.
Then, a remaining phase ramp in the object can optionally be computed by: - computing the sum of the calculated intensity along all frames - computing the shift of the center of mass of that intensity w/r to the array center - the phase ramp parameters are stored in the ptycho object, but not applied as the calculated
- Parameters:
obj – if True, after correcting both object and probe from the probe phase ramp, the object phase ramp is evaluated from the center of the calculated diffraction.
- pynx.ptycho.cl_operator.get_icalc(p: Ptycho, i, shift=False)#
Get the calculated intensity for frame i and a given Ptycho object. Note that this will reset the current psi, so should not be used during algorithms like DM which rely on the previous Psi value.
- Parameters:
p – the Ptycho object
i – the index of the desired frame
shift – if True, the array will be shifted so that the array is centered on the frame center (and not the corner).
- Returns:
the numpy array of the calculated intensity
Ptychography Runner class#
- class pynx.ptycho.runner.runner.PtychoRunner(argv, params, *args, **kwargs)#
Class to process a series of scans with a series of algorithms, given from the command-line.
- Parameters:
argv – the command-line parameters.
params – parameters for the optimization, with some default values.
args – ignored. For backward compatibility only
kwargs – ignored. For backward compatibility only
- check_params()#
Check if self.params includes a minimal set of valid parameters
Returns: Nothing. Will raise an exception if necessary
- check_params_beamline()#
Check if self.params includes a minimal set of valid parameters, specific to a beamline. Derived implementations can also set default values when appropriate.
Returns: Nothing. Will raise an exception if necessary
- classmethod make_parser(default_par, script_name='pynx-ptycho', description=None, epilog=None)#
Create the parser for the command-line arguments. This should be updated in derived objects.
- Parameters:
default_par – dictionary of default parameters
script_name – name of the script, e.g. ‘pynx-ptycho-cxi’. By default, will be read from the command-line arguments.
description – initial description string (top of help text)
epilog – epilog for the parser (usually with examples)
- Returns:
the parser
- parse_arg()#
Parses the arguments given on a command line
Returns: nothing
- parse_parameters_file(filename)#
Read parameters from a text file, written with one parameter per line. The parameters are converted to command-line arguments (e.g.
data = myfile.h5
is converted to--data myfile.h5
and appended to an argv-style list.- Parameters:
filename – the file to read the parameters from.
- Returns:
an argv-style list with the parameters parsed from the file
- print(*args, **kwargs)#
MPI-aware print function. Non-master processes will be muted :param args: args passed to print :param kwargs: kwrags passed to print :return: nothing
- process_scans()#
Run all the analysis on the supplied scan list, unless ‘help’ is given as a command-line argument.
- Returns:
Nothing
- exception pynx.ptycho.runner.runner.PtychoRunnerException#
- class pynx.ptycho.runner.runner.PtychoRunnerScan(params, scan, mpi_comm=None, timings=None)#
Abstract class to handle ptychographic data. Must be derived to be used. Only the load_scan() and load_data() functions need be derived.
- center_crop_data()#
Once the data has been loaded in self.load_data() [overloaded in child classes), this function can be called at the end of load_data to take care of centering the data (finding the center of diffraction) and cropping it with a size suitable for clFFT. Rebin is performed if necessary.
- Returns:
Nothing. self.iobs and self.dsize are updated. self.raw_data holds the raw data if needed for CXI export
- get_scan_prefix(run, scan=None)#
Get the save prefix given the parameters, scan and run numbers
- Parameters:
run – the run number
scan – the scan number - if None, self.scan is used
- init_run_number(_scan=None)#
Determine the current run number, depending on existing files :return: run0
- load_data()#
Loads data, using beamline-specific parameters. Abstract function, must be derived
Returns:
- load_data_post_process()#
Applies some post-processing to the input data, according to parameters. Also loads the mask. User-supplied mask is loaded if necessary.
This must be called at the end of load_data() :return:
- load_scan()#
Loads scan positions, using beamline-specific parameters. Abstract function, must be derived. This also filters the set of scan positions according to parameters (xyrange, monitor values, maxframe,…)
If MPI is used, only the master reads the scan positions. This automatically calls mpi_scan_split() at the end to split the scan if necessary (when MPI is used and mpi=’splitscan’ is used).
- Returns: Nothing. The scan positions, and the scan position indices ()
to be loaded for that runner are stored in self.x, self.y, self.imgn. If MPI is used, only the master should call this, and then call mpi_scan_split(). This is handled in Runner:process_scan().
- load_scan_post()#
Take care of frames with the direct beam illumination, mask positions. Should be called right after load_scan
- mpi_scan_split()#
This function is called after load_scan(), and will split the scan among all the MPI process. If MPI is not used, it does nothing and just passes the (x, y, imgn) values.
- Returns:
nothing. The x, y, imgn attributes are updated if necessary after splitting among MPI processes
- prepare()#
Prepare object and probe.
Returns: nothing. Adds self.obj0 and self.probe0
- prepare_processing_unit()#
Prepare processing unit (CUDA, OpenCL, or CPU). This must be called after load_scan so that the size of the dataset is known
Returns: nothing. Creates self.processing_unit, and adapts the stack size for CUDA
- print(*args, **kwargs)#
MPI-aware print function. Non-master processes will be muted :param args: args passed to print :param kwargs: kwrags passed to print :return: nothing
- print_probe_fwhm()#
Analyse probe shape and print estimated FWHM. Ignored for near field ptycho.
- Returns:
Nothing
- run(reuse_ptycho=False, share_probe_scan_data=None)#
Main work function, will run according to the set of algorithms specified in self.params.
- Parameters:
reuse_ptycho – if True, will reuse the previous Ptycho and PtychoData objects and skip some initialisation steps
share_probe_scan_data – used when processing jointly multiple scans with a shared probe
- Returns:
- run_algorithm(algo_string)#
Run a single or suite of algorithms in a given run.
- Parameters:
algo_string – a single or suite of algorithm steps to use, e.g. ‘ML**100’ or ‘analysis,ML**100,DM**100,nbprobe=3,DM**100’ or ‘analysis,ML**100*DM**100,nbprobe=3,DM**100’
- Returns:
Nothing
- save(run, stepnum=None, algostring=None, scan=None)#
Save the result of the optimization, and (if self.params[‘saveplot’] is True) the corresponding plot. This is an internal function.
- Parameters:
run – the run number (integer)
stepnum – the step number in the set of algorithm steps
algostring – the string corresponding to all the algorithms ran so far, e.g. ‘100DM,100AP,100ML’
scan – the scan number - normally self.scan, but can be used when processing simultaneously multiple scans
- Returns:
- save_data_cxi(crop=True, verbose=False, **kwargs)#
Save the scan data using the CXI format (see http://cxidb.org) :param crop: if True, only the already-cropped data will be save. Otherwise, the original raw data is saved.
Returns:
- save_movie()#
Create a movie of the scan with all scan positions and diffraction frames. Requires matplotlib and ffmpeg. :return:
- save_plot(run, stepnum=None, algostring=None, display_plot=False)#
Save the plot to a png file.
- Parameters:
run – the run number (integer)
stepnum – the step number in the set of algorithm steps
algostring – the string corresponding to all the algorithms ran so far, e.g. ‘100DM,100AP,100ML’
display_plot – if True, the saved plot will also be displayed
- Returns:
- test_no_rerun(verbose=True)#
Test if this scan # should be skipped; :param scan: the scan number :return: True if the scan should be skipped, False otherwise
Ptychography Analysis#
- pynx.ptycho.analysis.modes(d, pixel_size, do_plot=True, show_plot=True, verbose=False)#
Determine complex modes for a given object or probe.
- Parameters:
d – the probe or object to analyse, with the modes along the first axis/
pixel_size – the pixel size in meters
do_plot – plot if True (default)
show_plot – if True, the plot will be done in a visible figure, otherwise in an offline one. Ignored unless do_plot is True.
- Returns:
a tuple of (orthogonalized modes with the same shape as the input, a vector of the relative intensities, figure)
- pynx.ptycho.analysis.probe_fwhm(probe, pixel_size, verbose=True)#
Analyse probe shape and estimated FWHM using different methods:
Full width at half maximum
Full width at 20%
Full width at half maximum using a statistical gaussian analysis - the width is the same as a gaussian with a maximum intensity of 1 and the same integrated intensity.
- Parameters:
probe – the probe to analyse. If number of dimensions is 3, only the first mdoe is analysed
pixel_size – the pixel size in meters
verbose – if True (the default), will print the corresponding sizes
- Returns:
((fwhm_x, fwhm_y), (fwhm_x@20%, fwhm_y@20%), (fwhm_gauss_stat, fwhm_gauss_stat_x, fwhm_gauss_stat_y))
- pynx.ptycho.analysis.probe_propagate(probe, propagation_range, pixel_size, wavelength, do_plot=True, show_plot=True, fig_size=(15, 10))#
Propagate a probe along a given range, plot it and return the probe at the focus (position where the size is minimal)
- Parameters:
probe – the 2D complex probe to propagate
propagation_range – either given as a range/array of dz values (in meters), or a tuple with (zmin,zmax,nb)
pixel_size – the probe pixel size, in meters
do_plot – plot if True (default)
show_plot – if True, the plot will be done in a visible figure, otherwise in an offline one. Ignored unless do_plot is True
fig_size – the figure size, default = (15,10)
- Returns:
the 3d array of the propagated probe along z, with size (nz, ny, nx)
the z coordinates along the 3d array
the index of the found focus point along the z direction
the figure if it was plotted
- Return type:
A tuple with
Examples#
Operator-based API, far field#
########################################################################
#
# Example of the ptychograpic reconstruction using OpenCL on simulated data
# (c) ESRF 2017-present
# Authors: Vincent Favre-Nicolin <favre@esrf.fr>
#
########################################################################
import timeit
from pylab import *
if False:
from pynx.processing_unit import default_processing_unit
# This can be used to select the GPU and/or the language, and must be called before other pynx imports
# Otherwise a default GPU will be used according to its speed, which is usually sufficient
default_processing_unit.select_gpu(language='OpenCL', gpu_name='R9')
# The following statement will import either CUDA or OpenCL operators
from pynx.ptycho import *
from pynx.ptycho import simulation
##################
# Simulation of the ptychographic data:
n = 256
pixel_size_detector = 55e-6
wavelength = 1.5e-10
detector_distance = 1
obj_info = {'type': 'phase_ampl', 'phase_stretch': pi / 2, 'alpha_win': .2}
# probe_info = {'type': 'focus', 'aperture': (30e-6, 30e-6), 'focal_length': .08, 'defocus': 100e-6, 'shape': (n, n)}
probe_info = {'type': 'gauss', 'sigma_pix': (40, 40), 'defocus': 100e-6, 'shape': (n, n)}
# 50 scan positions correspond to 4 turns, 78 to 5 turns, 113 to 6 turns
scan_info = {'type': 'spiral', 'scan_step_pix': 30, 'n_scans': 120}
data_info = {'num_phot_max': 1e9, 'bg': 0, 'wavelength': wavelength, 'detector_distance': detector_distance,
'detector_pixel_size': pixel_size_detector, 'noise': 'poisson'}
# Initialisation of the simulation with specified parameters
s = simulation.Simulation(obj_info=obj_info, probe_info=probe_info, scan_info=scan_info, data_info=data_info)
s.make_data()
# Positions from simulation are given in pixels
posx, posy = s.scan.values
ampl = s.amplitude.values # square root of the measured diffraction pattern intensity
pixel_size_object = wavelength * detector_distance / pixel_size_detector / n
##################
# Size of the reconstructed object (obj)
nyo, nxo = shape.calc_obj_shape(posx, posy, ampl.shape[1:])
# Initial object
# obj_init_info = {'type':'flat','shape':(nx,ny)}
obj_init_info = {'type': 'random', 'range': (0, 1, 0, 0.5), 'shape': (nyo, nxo)}
# Initial probe
probe_init_info = {'type': 'focus', 'aperture': (20e-6, 20e-6), 'focal_length': .08, 'defocus': 50e-6, 'shape': (n, n)}
data_info = {'wavelength': wavelength, 'detector_distance': detector_distance,
'detector_pixel_size': pixel_size_detector}
init = simulation.Simulation(obj_info=obj_init_info, probe_info=probe_init_info, data_info=data_info)
init.make_obj()
init.make_probe()
data = PtychoData(iobs=ampl ** 2, positions=(posx * pixel_size_object, posy * pixel_size_object), detector_distance=1,
mask=None, pixel_size_detector=55e-6, wavelength=1.5e-10)
p = Ptycho(probe=s.probe.values, obj=init.obj.values, data=data, background=None)
# Initial scaling is important to avoid overflows during ML
p = ScaleObjProbe(verbose=True) * p
p = DM(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 100 * p
p = ML(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 40 * p
if True:
print("###########################################################################################################")
# Add probe modes
pr = p.get_probe()
ny, nx = pr.shape[-2:]
nb_probe = 3
pr1 = np.empty((nb_probe, ny, nx), dtype=np.complex64)
pr1[0] = pr[0]
for i in range(1, nb_probe):
n = abs(pr).mean()
pr1[i] = np.random.uniform(0, n, (ny, nx)) * exp(1j * np.random.uniform(0, 2 * np.pi, (ny, nx)))
p.set_probe(pr1)
p = DM(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 100 * p
p = AP(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 40 * p
p = ML(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 100 * p
if False:
# Timing vs stack size
n = 50
vx = []
vy = []
for stack_size in range(1, 32 + 1, 1):
default_processing_unit.set_stack_size(stack_size)
p = DM(update_object=True, update_probe=True) ** 10 * p
t0 = timeit.default_timer()
p = DM(update_object=True, update_probe=True) ** n * p
dt = (timeit.default_timer() - t0) / n
print("DM dt/cycle=%5.3fs [stack_size=%2d]" % (dt, stack_size))
vx.append(stack_size)
vy.append(dt)
figure()
plot(vx, vy, '-')
xlabel('stack size')
ylabel('Time for a DM cycle (s)')
if False:
# Export data and object, probe to CXI files
p.save_obj_probe_cxi('obj_probe.cxi')
save_ptycho_data_cxi('data.cxi', ampl ** 2, pixel_size_detector, wavelength, detector_distance,
posx * pixel_size_object, posy * pixel_size_object, z=None, monitor=None,
mask=None, instrument='simulation', overwrite=True)
Operator-based API, near field#
########################################################################
#
# Example of the ptychograpic reconstruction using OpenCL on simulated data
# (c) ESRF 2017-present
# Authors: Vincent Favre-Nicolin <favre@esrf.fr>
#
########################################################################
import timeit
from pylab import *
from pynx.ptycho import simulation, shape
# from pynx.ptycho import *
from pynx.ptycho.ptycho import *
from pynx.ptycho.cpu_operator import * # Use only CPU operators ?
from pynx.ptycho.operator import * # Use CUDA > OpenCL > CPU operators, as available
import pynx.ptycho.cpu_operator as cpuop
import pynx.ptycho.cl_operator as clop
##################
detector_distance = 1.5
wavelength = 1.5e-10
pixel_size_detector = 1e-6
# Simulation of the ptychographic data:
n = 256
obj_info = {'type': 'phase_ampl', 'phase_stretch': pi / 2, 'alpha_win': .2}
probe_info = {'type': 'near_field', 'aperture': (80e-6, 80e-6), 'defocus': 0.5, 'shape': (n, n)}
# 50 scan positions correspond to 4 turns, 78 to 5 turns, 113 to 6 turns
scan_info = {'type': 'spiral', 'scan_step_pix': 20, 'n_scans': 120}
data_info = {'num_phot_max': 1e9, 'bg': 0, 'wavelength': wavelength, 'detector_distance': detector_distance,
'detector_pixel_size': pixel_size_detector, 'noise': 'poisson', 'near_field': True}
# Initialisation of the simulation
s = simulation.Simulation(obj_info=obj_info, probe_info=probe_info, scan_info=scan_info, data_info=data_info)
s.make_data()
# Positions from simulation are given in pixels
posx, posy = s.scan.values
ampl = s.amplitude.values # square root of the measured diffraction pattern intensity
##################
# Size of the reconstructed object (obj)
nyo, nxo = shape.calc_obj_shape(posx, posy, ampl.shape[1:])
# Initial object
# obj_init_info = {'type':'flat','shape':(nx,ny)}
obj_init_info = {'type': 'random', 'range': (0, 1, 0, 0.5), 'shape': (nyo, nxo)}
# Initial probe
probe_init_info = {'type': 'near_field', 'aperture': (90e-6, 90e-6), 'defocus': 0.3, 'shape': (n, n)}
init = simulation.Simulation(obj_info=obj_init_info, probe_info=probe_init_info, data_info=data_info)
init.make_obj()
init.make_probe()
data = PtychoData(iobs=ampl ** 2, positions=(posx * pixel_size_detector, posy * pixel_size_detector),
detector_distance=detector_distance, mask=None,
pixel_size_detector=pixel_size_detector, wavelength=wavelength, near_field=True)
p = Ptycho(probe=init.probe.values, obj=init.obj.values, data=data, background=None)
# Initial scaling is important to avoid overflows during ML
p = ScaleObjProbe(verbose=True) * p
p = DM(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 100 * p
p = AP(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 40 * p
p = ML(update_object=True, update_probe=True, calc_llk=20, show_obj_probe=20) ** 40 * p
if False:
# Timing vs stack size
n = 50
vx = []
vy = []
for stack_size in range(1, 32 + 1, 1):
default_processing_unit.set_stack_size(stack_size)
p = DM(update_object=True, update_probe=True) ** 10 * p
t0 = timeit.default_timer()
p = DM(update_object=True, update_probe=True) ** n * p
dt = (timeit.default_timer() - t0) / n
print("DM dt/cycle=%5.3fs [stack_size=%2d]" % (dt, stack_size))
vx.append(stack_size)
vy.append(dt)
figure()
plot(vx, vy, '-')
xlabel('stack size')
ylabel('Time for a DM cycle (s)')