pynx.cdi: Coherent Diffraction Imaging

Description

This modules provides algorithm for 2D and 3D reconstruction of single objects using several algorithms:

  • hybrid input-output (HIO)

  • error reduction (ER)

  • relaxed averaged alternating reflections (RAAR),

  • maximum likelihood conjugate gradient minimization…

  • partial coherence

  • etc…

The calculations use an ‘operator’ approach, where operations on a given cdi dataset can be simply written by multiplying it by the corresponding operator (FT, projection,..) or by a series of operators.

API Reference

Note that the Python API uses an ‘operator’ approach, to enable writing complex operations in a more mathematical and natural way.

CDI base classes

This is the CDI base classes, which can be used with operators

class pynx.cdi.cdi.CDI(iobs, support=None, obj=None, mask=None, pixel_size_detector=None, wavelength=None, detector_distance=None, crop=None, bin=None)

Reconstruction class for two or three-dimensional coherent diffraction imaging data.

Constructor. All arrays are assumed to be centered in (0,0) to avoid fftshift

Parameters:
  • iobs

    2D/3D observed diffraction data (intensity). Assumed to be corrected and following Poisson statistics, will be converted to float32. Dimensions should be divisible by 4 and have a prime factor decomposition up to 7. Internally, the following special values are used:

    • values<=-1e19 are masked. Among those, values in ]-1e38;-1e19] are estimated values, stored as -(iobs_est+1)*1e19, which can be used to make a loose amplitude projection. Values <=-1e38 are masked (no amplitude projection applied), just below the minimum float32 value

    • -1e19 < values <= 1 are observed but used as free pixels If the mask is not supplied, then it is assumed that the above special values are used.

  • support – initial support in real space (1 = inside support, 0 = outside)

  • obj – initial object. If None, it should be initialised later.

  • mask – mask for the diffraction data (0: valid pixel, >0: masked)

  • pixel_size_detector – detector pixel size (meters)

  • wavelength – experiment wavelength (meters)

  • detector_distance – detector distance (meters)

copy(copy_history_llk=False, copy_free_pixels=False)

Creates a copy (without any reference passing) of this object, unless copy_obj is False.

Parameters:
  • copy_history_llk – if True, the history, cycle number, llk, etc.. are copied to the new object.

  • copy_free_pixels – if True, copy the distribution of free pixels to the new CDI object. Otherwise the new CDI object does not have free pixels initialised.

Returns:

a copy of the object.

get_bin(dim3=False)

Get the binning parameter as an int32 array, one value per dimension.

Parameters:

dim3 – this can b used to force the number of elements to 3, for convenience in GPU kernels.

Returns:

the bin parameter

get_crop(dim3=False)

Get the cropping parameter as an int32 array, one value per dimension.

Parameters:

dim3 – this can b used to force the number of elements to 3, for convenience in GPU kernels.

Returns:

the bin parameter

get_iobs(shift=False)

Get the observed intensity data array.

Parameters:

shift – if True, the data array will be fft-shifted so that the center of the data is in the center of the array, rather than in the corner (the default).

Returns:

the 2D or 3D CDI numpy data array

get_llk(noise=None, normalized=True)

Get the normalized log-likelihoods, which should converge to 1 for a statistically ideal fit.

Parameters:
  • noise – either ‘gaussian’, ‘poisson’ or ‘euclidian’, will return the corresponding log-likelihood.

  • normalized – if True, will return normalized values so that the llk from a statistically ideal model should converge to 1

Returns:

the log-likelihood, or if noise=None, a tuple of the three (poisson, gaussian, euclidian) log-likelihoods.

get_llkn()

Get the poisson normalized log-likelihood, which should converge to 1 for a statistically ideal fit to a Poisson-noise dataset.

Returns:

the normalized log-likelihood, for Poisson noise.

get_obj(shift=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:

shift – if True, the data array will be fft-shifted so that the center of the data is in the center of the array, rather than in the corner (the default).

Returns:

the 2D or 3D CDI numpy data array

get_support(shift=False)

Get the support array. This will automatically get the latest data, either from GPU or from the host memory, depending where the last changes were made.

Parameters:

shift – if True, the data array will be fft-shifted so that the center of the data is in the center of the array, rather than in the corner (the default).

Returns:

the 2D or 3D CDI numpy data array

get_upsample(dim3=False)

Get the binning parameter as an int32 array, one value per dimension.

Parameters:

dim3 – this can b used to force the number of elements to 3, for convenience in GPU kernels.

Returns:

the bin parameter

get_x_y_z()

Get 1D arrays of x and y (z if 3d) coordinates, taking into account the pixel size. The arrays are centered at (0,0) - i.e. with the origin in the corner for FFT puroposes. x is an horizontal vector and y vertical, and (if 3d) z along third dimension.

Returns:

a tuple (x, y) or (x, y, z) of 1D numpy arrays

in_object_space()
Returns:

True if the current obj array is in object space, False otherwise.

init_free_pixels(ratio=0.05, island_radius=3, exclude_zone_center=0.05, verbose=False, coords=None)

Random selection of ‘free’ pixels in the observed intensities, to be used as unbiased indicator. The free pixel iobs values are modified as iobs_free = -iobs - 1

Parameters:
  • ratio – (approximate) relative number of pixels to be included in the free mask

  • island_radius – free island radius, to avoid pixel correlation due to finite object size

  • exclude_zone_center – the relative radius of the zone to be excluded near the center

  • verbose – if True, be verbose

  • coords – instead of generating random coordinates, these can be given as a tuple of (ix, iy[, iz]). All coordinates should be at least island_radius far from the borders, and these coordinates are relative to the array centre, i.e. within [-size/2 + island_radius ; size/2 - island_radius] along each axis.

Returns:

a tuple (ix, iy[, iz]) of the pixel coordinates of the islands, before taking into account any mask.

init_psf(model='pseudo-voigt', fwhm=1, eta=0.1, psf_f=None)

Initialise the point-spread function to model the partial coherence, using either a Lorentzian, Gaussian or pseudo-Voigt function

Parameters:
  • model – “lorentzian”, “gaussian” or “pseudo-voigt”, or None to deactivate

  • fwhm – the full-width at half maximum, in pixels

  • eta – the eta parameter for the pseudo-voigt

  • psf_f – this can be used to supply an array for the PSF. This should be the half-Hermition which is the result of the real2complex FT of the detector-space PSF. It should be centered at the origin (fft-shifted). If this is given all other arguments are ignored.

Returns:

nothing

reset_history()

Reset history, and set current cycle to zero :return: nothing

save_data_cxi(filename, sample_name=None, experiment_id=None, instrument=None, process_parameters=None)

Save the diffraction data (observed intensity, mask) to an HDF% CXI file. :param filename: the file name to save the data to :param sample_name: optional, the sample name :param experiment_id: the string identifying the experiment, e.g.: ‘HC1234: Siemens star calibration tests’ :param instrument: the string identifying the instrument, e.g.: ‘ESRF id10’ :param process_parameters: a dictionary of parameters which will be saved as a NXcollection :return: Nothing. a CXI file is created

save_obj_cxi(filename, sample_name=None, experiment_id=None, instrument=None, note=None, crop=0, save_psf=True, process_notes=None, process_parameters=None, append=False, shift_phase_zero=True, **kwargs)

Save the result of the optimisation (object, support) to an HDF5 CXI file. :param filename: the file name to save the data to :param sample_name: optional, the sample name :param experiment_id: the string identifying the experiment, e.g.: ‘HC1234: Siemens star calibration tests’ :param instrument: the string identifying the instrument, e.g.: ‘ESRF id10’ :param note: a string with a text note giving some additional information about the data, a publication… :param crop: integer, if >0, the object will be cropped arount the support, plus a margin of ‘crop’ pixels. :param save_psf: if True, also save the psf (if present) :param process_notes: a dictionary which will be saved in ‘/entry_N/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}} :param process_parameters: a dictionary of parameters which will be saved as a NXcollection :param append: if True and the file already exists, the result will be saved in a new entry. :param shift_phase_zero: if True, try to center the object phase around zero inside the support :return: Nothing. a CXI file is created

set_bin(bin_f, obj_support=True, update_nb_points=True)

Change the binning parameters, either as a single integer (same value for all dimensions) or as a tuple/list of values. The observed intensity array size is then binned (and the original array is kept), which results in an object with the same resolution, but a smaller real-space extension, i.e. a lower oversampling (so support update should be done more carefully). When binning the Iobs array, the pixels tagged for free LLK calculation are handled. The only values used should be 1 or 2.

This automatically disables crop and upsample parameters.

Parameters:
  • bin_f – if bin_f is not None or 1, the iobs array size is binned by this factor, keeping the object pixel size constant but with a smaller extent. A list/tuple of integers can also be supplied to use a different binning for each dimension.

  • obj_support – if True, the object and support arrays will be modified according to the bin change in the iobs array. Using obj_support=False allows to perform separately the object and support scaling using a GPU operator.

Returns:

nothing, the binning parameter is stored, the iobs array is binned, the object array is resized.

set_crop(crop_f, obj_support=True, update_nb_points=True)

Change the cropping parameters, either as a single integer (same value for all dimensions) or as a tuple/list of values. The observed intensity array size is then cropped (and the original array is kept), which results in an object with a lower resolution, but the same extent. The only values used should be 1 or 2.

This automatically disables bin and upsample parameters.

Parameters:
  • crop_f – if crop_f is not None or 1, the iobs array is cropped by this factor, also dividing the object size by the same amount, with the object pixel size multiplied by the factor. A list/tuple of integers can also be supplied to use a different binning for each dimension.

  • obj_support – if True, the object and support arrays will be modified according to the bin change in the iobs array. Using obj_support=False allows to perform separately the object and support scaling using a GPU operator.

Returns:

nothing, the crop parameter is stored and will be used for computations.

set_free_pixel_mask(m, verbose=False, shift=False)

Set the free pixel mask. Assumes the free pixel mask is correct (excluding center, bad pixels)

Parameters:
  • m – the boolean mask (1=masked, 0 otherwise)

  • shift – if True, the mask array will be fft-shifted so that the center of the stored data is in the corner of the array. [default: the array is already shifted]

Returns:

nothing

set_iobs(iobs, shift=False)

Set the observed intensity data array.

Parameters:
  • iobs – the 2D or 3D CDI numpy data array (float32 numpy array)

  • shift – if True, the data array will be fft-shifted so that the center of the stored data is in the corner of the array. [default: the array is already shifted]

Returns:

nothing

set_mask(mask, shift=False)

Set the mask data array. Note that since the mask is stored by setting observed intensities of masked pixels to negative values, it is not possible to un-mask pixels.

Parameters:
  • obj – the 2D or 3D CDI mask array

  • shift – if True, the data array will be fft-shifted so that the center of the stored data is in the corner of the array. [default: the array is already shifted]

Returns:

nothing

set_obj(obj, shift=False)

Set the object data array. Assumed to be in object (not Fourier) space

Parameters:
  • obj – the 2D or 3D CDI numpy data array (complex64 numpy array)

  • shift – if True, the data array will be fft-shifted so that the center of the stored data is in the corner of the array. [default: the array is already shifted]

Returns:

nothing

set_support(support, shift=False)

Set the support data array.

Parameters:
  • obj – the 2D or 3D CDI numpy data array (complex64 numpy array)

  • shift – if True, the data array will be fft-shifted so that the center of the stored data is in the corner of the array. [default: the array is already shifted]

Returns:

nothing

set_upsample(upsample_f)

Change the upsample parameters, either as a single integer (same value for all dimensions) or as a tuple/list of values. The observed intensity array size is not affected by this parameter, but the extent of the object will be multiplied by this factor. The upsampled calculated intensities will be binned for comparison with the observed intensity array. This can be used when the oversampling ratio is not large enough for a reconstruction [Maddali et al, Phys. Rev. A. 99 (2019), 053838]

This does not alter the observed intensity array, or the object array which is accessible by get_obj(), which will keep the dimensions corresponding to the original Iobs array.

This automatically disables bin and crop parameters.

Parameters:

upsample_f – if upsample_f is not None or 1, the object array size is internally extended by this factor, keeping the object pixel size constant. A list/tuple of integers can also be supplied to use a different upsampling for each dimension. This does not affect the resolution of the reconstructed object, but simulates an increase of the oversampling ratio by the upsampling factor.

Returns:

nothing, the upsampling parameter is stored and will be used for computations.

update_history(mode='algorithm', verbose=False, **kwargs)

Update the history record. :param mode: either ‘llk’ (will record new log-likelihood, nb_photons_calc, average value..), or ‘algorithm’ (will only update the algorithm) or ‘support’. :param verbose: if True, print some info about current process :param kwargs: other parameters to be recorded, e.g. support_threshold=xx, dt= :return: nothing

exception pynx.cdi.cdi.CDIOperatorException
class pynx.cdi.cdi.OperatorCDI(lazy=False)

Base class for an operator on CDI objects, not requiring a processing unit.

timestamp_increment(cdi)

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:

exception pynx.cdi.cdi.SupportTooLarge
exception pynx.cdi.cdi.SupportTooSmall
pynx.cdi.cdi.calc_throughput(p: CDI | None = None, cxi=None, verbose=False)

Analyse the throughput after a series of algorithms, either from a CDI object or from a CXI file. This does not take into account operations for support update, or partial coherence (PSF). :param p: the CDI object the timings are extracted from. :param cxi: the CXI file the history of cycles will be obtained from :param verbose: if True, print average throughput per algorithm steps :return: the average throughput in Gbyte/s

pynx.cdi.cdi.save_cdi_data_cxi(filename, iobs, wavelength=None, detector_distance=None, pixel_size_detector=None, mask=None, sample_name=None, experiment_id=None, instrument=None, note=None, iobs_is_fft_shifted=False, process_parameters=None)

Save the diffraction data (observed intensity, mask) to an HDF5 CXI file, NeXus-compatible. :param filename: the file name to save the data to :param iobs: the observed intensity :param wavelength: the wavelength of the experiment (in meters) :param detector_distance: the detector distance (in meters) :param pixel_size_detector: the pixel size of the detector (in meters) :param mask: the mask indicating valid (=0) and bad pixels (>0) :param sample_name: optional, the sample name :param experiment_id: the string identifying the experiment, e.g.: ‘HC1234: Siemens star calibration tests’ :param instrument: the string identifying the instrument, e.g.: ‘ESRF id10’ :param iobs_is_fft_shifted: if true, input iobs (and mask if any) have their origin in (0,0[,0]) and will be shifted back to centered-versions before being saved. :param process_parameters: a dictionary of parameters which will be saved as a NXcollection :return: Nothing. a CXI file is created

CDI Operators

This section lists the operators, which can be imported automatically using from pynx.cdi 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.cdi.cl_operator.ApplyAmplitude(calc_llk=False, zero_mask=False, scale_in=1, scale_out=1, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=False, psf_filter=None)

Apply the magnitude from an observed intensity, keep the phase.

Parameters:
  • calc_llk – if true, the log-likelihood will be calculated and stored in the object

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • 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

  • confidence_interval_factor – a relaxation factor, with the projection of calculated amplitude being done towards the limit of the poisson confidence interval. A value of 1 corresponds to a 50% confidence interval, a value of 0 corresponds to a strict observed amplitude projection. [EXPERIMENTAL]

  • confidence_interval_factor_mask_max (confidence_interval_factor_mask_min,) – For masked pixels where a value has been estimated (e.g. with InterpIobsMask()), a confidence interval can be given as a factor to be applied to the interpolated observed intensity. This corresponds to values stored between -1e19 and -1e38. [EXPERIMENTAL]

  • update_psf – if True, will update the PSF convolution kernel using the Richard-Lucy deconvolution approach. If there is no PSF, it will be automatically initialised.

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.AutoCorrelationSupport(threshold=0.2, verbose=False, lazy=False, scale=False)

Operator to calculate an initial support from the auto-correlation function of the observed intensity. The object will be multiplied by the resulting support.

Parameters:
  • threshold – pixels above the autocorrelation maximum multiplied by the threshold will be included in the support. This can either be a float, or a range tuple (min,max) between which the threshold value will be randomly chosen every time the operator is applied.

  • verbose – if True, print info about the result of the auto-correlation

  • lazy – if True, this will be queued for later execution in the cdi object

  • scale – if True, or a string is given, will also apply ScaleObj(method=scale)

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.CF(positivity=False, calc_llk=False, nb_cycle=1, show_cdi=False, fig_num=-1, zero_mask=False, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=0, psf_filter=None, plot_axis=0)

Charge flipping cycle

Parameters:
  • positivity – apply a positivity restraint

  • calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle

  • nb_cycle – the number of cycles to perform

  • show_cdi – 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 intensity, as for ShowCDI()

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • confidence_interval_factor – a relaxation factor, with the projection of calculated amplitude being done towards the limit of the poisson confidence interval. A value of 1 corresponds to a 50% confidence interval, a value of 0 corresponds to a strict observed amplitude projection. [EXPERIMENTAL]

  • confidence_interval_factor_mask_max (confidence_interval_factor_mask_min,) – For masked pixels where a value has been estimated (e.g. with InterpIobsMask()), a confidence interval can be given as a factor to be applied to the interpolated observed intensity. This corresponds to values stored between -1e19 and -1e38. [EXPERIMENTAL]

  • update_psf – if >0, will update the partial coherence psf every update_psf cycles.

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.DetwinHIO(detwin_axis=0, nb_cycle=10, beta=0.9, positivity=False, zero_mask=False)

HIO cycles with a temporary halved support

Constructor for the DetwinHIO operator

Parameters:
  • detwin_axis – axis along which the detwinning will be performed. If None, a random axis is chosen

  • nb_cycle – number of cycles to perform while using a halved support

  • beta – the beta value for the HIO operator

  • positivity – True or False

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.DetwinRAAR(detwin_axis=0, nb_cycle=10, beta=0.9, positivity=False, zero_mask=False)

RAAR cycles with a temporary halved support

Constructor for the DetwinRAAR operator

Parameters:
  • detwin_axis – axis along which the detwinning will be performed. If None, a random axis is chosen

  • nb_cycle – number of cycles to perform while using a halved support

  • beta – the beta value for the HIO operator

  • positivity – True or False

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.ER(positivity=False, calc_llk=False, nb_cycle=1, show_cdi=False, fig_num=-1, zero_mask=False, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=0, psf_filter=None, plot_axis=0)

Error reduction cycle

Parameters:
  • positivity – apply a positivity restraint

  • calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle

  • nb_cycle – the number of cycles to perform

  • show_cdi – 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 intensity, as for ShowCDI()

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • confidence_interval_factor – a relaxation factor, with the projection of calculated amplitude being done towards the limit of the poisson confidence interval. A value of 1 corresponds to a 50% confidence interval, a value of 0 corresponds to a strict observed amplitude projection. [EXPERIMENTAL]

  • confidence_interval_factor_mask_max (confidence_interval_factor_mask_min,) – For masked pixels where a value has been estimated (e.g. with InterpIobsMask()), a confidence interval can be given as a factor to be applied to the interpolated observed intensity. This corresponds to values stored between -1e19 and -1e38. [EXPERIMENTAL]

  • update_psf – if >0, will update the partial coherence psf every update_psf cycles.

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.EstimatePSF(*args, **kwargs)

Estimate the Point Spread Function. [OBSOLETE, replaced by InitPSF]

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.FT(scale=True)

Forward Fourier transform.

Parameters:

scale – if True, the Fourier transform will be normalised, so that the transformed array L2 norm will remain constant (by dividing the output by the square root of the object’s size). If False or None, the array norm will not be changed. If a scalar is given, the output array is multiplied by it.

op(cdi)

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).

class pynx.cdi.cl_operator.FourierApplyAmplitude(calc_llk=False, zero_mask=False, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=False, psf_filter=None, obj_stats=False)

Fourier magnitude operator, performing a Fourier transform, the magnitude projection, and a backward FT.

Parameters:
  • calc_llk – if True, the log-likelihood will be calculated while in diffraction space.

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • obj_stats – if True, will call ObjSupportStats at the end

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

op(cdi)

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).

class pynx.cdi.cl_operator.FreePU(processing_unit=None, lazy=False)

Operator freeing OpenCL memory. The FFT plan/app in self.processing_unit is removed, as well as any OpenCL cla.Array attribute in the supplied CDI object.

The latest object and support data is retrieved from GPU memory

op(cdi)

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(cdi)

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.cdi.cl_operator.GPS(inertia=0.05, t=1.0, s=0.9, sigma_f=0, sigma_o=0, positivity=False, calc_llk=False, nb_cycle=1, show_cdi=False, fig_num=-1, zero_mask=False, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=0, psf_filter=None, plot_axis=0)

GPS cycle, according to Pham et al [2019]

Parameters:
  • inertia – inertia parameter (sigma in original Pham2019 article)

  • t – t parameter

  • s – s parameter

  • sigma_f – Fourier-space smoothing kernel width, in Fourier-space pixel units

  • sigma_o – object-space smoothing kernel width, in object-space pixel units

  • positivity – if True, apply a positivity restraint

  • calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle

  • nb_cycle – the number of cycles to perform

  • show_cdi – 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 intensity, as for ShowCDI()

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • confidence_interval_factor – a relaxation factor, with the projection of calculated amplitude being done towards the limit of the poisson confidence interval. A value of 1 corresponds to a 50% confidence interval, a value of 0 corresponds to a strict observed amplitude projection. [EXPERIMENTAL]

  • confidence_interval_factor_mask_max (confidence_interval_factor_mask_min,) – For masked pixels where a value has been estimated (e.g. with InterpIobsMask()), a confidence interval can be given as a factor to be applied to the interpolated observed intensity. This corresponds to values stored between -1e19 and -1e38. [EXPERIMENTAL]

  • update_psf – if >0, will update the partial coherence psf every update_psf cycles.

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.HIO(beta=0.9, positivity=False, calc_llk=False, nb_cycle=1, show_cdi=False, fig_num=-1, zero_mask=False, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=0, psf_filter=None, plot_axis=0)

Hybrid Input-Output reduction cycle

Parameters:
  • positivity – apply a positivity restraint

  • calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle

  • nb_cycle – the number of cycles to perform

  • show_cdi – 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 intensity, as for ShowCDI()

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • confidence_interval_factor – a relaxation factor, with the projection of calculated amplitude being done towards the limit of the poisson confidence interval. A value of 1 corresponds to a 50% confidence interval, a value of 0 corresponds to a strict observed amplitude projection. [EXPERIMENTAL]

  • confidence_interval_factor_mask_max (confidence_interval_factor_mask_min,) – For masked pixels where a value has been estimated (e.g. with InterpIobsMask()), a confidence interval can be given as a factor to be applied to the interpolated observed intensity. This corresponds to values stored between -1e19 and -1e38. [EXPERIMENTAL]

  • update_psf – if >0, will update the partial coherence psf every update_psf cycles.

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.IFT(scale=True)

Inverse Fourier transform

Parameters:

scale – if True, the Fourier transform will be normalised, so that the transformed array L2 norm will remain constant (by multiplying the output by the square root of the object’s size). If False or None, the array norm will not be changed. If a scalar is given, the output array is multiplied by it.

op(cdi)

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).

class pynx.cdi.cl_operator.InitPSF(model='pseudo-voigt', fwhm=1, eta=0.05, psf=None, filter=None)

Initialise the point-spread function kernel to model partial coherence.

Initialise the point-spread function to model the partial coherence, using either a Lorentzian, Gaussian or Pseuo-Voigt function.

Parameters:
  • model – “lorentzian”, “gaussian” or “pseudo-voigt”. The default is a pseudo-Voigt, as it allows a relatively sharp peak while still keeping some tails which allow the psf kernel to be updated.

  • fwhm – the full-width at half maximum, in pixels

  • eta – eta value for the pseudo-Voigt (default 0.01)

  • psf – an array of the PSF can be supplied. In that case the other parameters are ignored. The psf can be smaller than the iobs array size, and will be resized, normalised and transformed to the reciprocal space kernel used internally. The psf should be centred on the array centre, and will be fft-shifted automatically.

  • filter – None, “hann” or “tukey” - filter for the initial PSF update. This is not used if the psf array is given as a parameter.

Returns:

nothing. This initialises cdi._cl_psf_f, and copies the array to cdi._psf_f

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.InterpIobsMask(d=8, n=4)

Interpolate masked pixels observed intensity using inverse distance weighting

Parameters:
  • d – the half-distance of the interpolation, which will be done for pixel i from i-d to i+d along each dimension

  • n – the weighting will be calculated as 1/d**n

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.LLK(scale=1.0)

Log-likelihood reduction kernel. This is a reduction operator - it will write llk as an argument in the cdi object. If it is applied to a CDI instance in object space, a FT() and IFT() will be applied to perform the calculation in diffraction space. This collects log-likelihood for Poisson, Gaussian and Euclidian noise models, and also computes the total calculated intensity (including in masked pixels).

Parameters:

scale – the scale factor to be applied to the calculated amplitude before evaluating the log-likelihood. The calculated amplitudes are left unmodified.

op(cdi)

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).

class pynx.cdi.cl_operator.LLKSupport(processing_unit=None, lazy=False)

Support log-likelihood reduction kernel. Can only be used when cdi instance is object space. This is a reduction operator - it will write llk_support as an argument in the cdi object, and return cdi.

op(cdi)

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).

class pynx.cdi.cl_operator.ML(reg_fac=0.01, nb_cycle=1, calc_llk=False, show_cdi=False, fig_num=-1, plot_axis=0)

Maximum likelihood conjugate gradient minimization

Parameters:
  • reg_fac

  • nb_cycle – the number of cycles to perform

  • 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_cdi – 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 intensity, as for ShowCDI()

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.ObjConvolve(sigma=1)

Gaussian convolution of the object, produces a new array with the convoluted amplitude of the object.

op(cdi)

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).

class pynx.cdi.cl_operator.PRTF(fig_num=None, file_name=None, nb_shell=None, fig_title=None)

Operator to compute the Phase Retrieval Transfer Function. When applied to a CDI object, it stores the result in it as cdi.prtf, cdi.prtf_freq, cdi.prtf_nyquist, cdi.prtf_nb

Parameters:
  • fig_num – the figure number to display the PRTF.

  • file_name – if given, the PRTF figure will be saved to this file (should end in .png or .pdf).

  • nb_shell – the number of shell in which to compute the PRTF. By default the shell thickness is 2 pixels

  • fig_title – the figure title

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.RAAR(beta=0.9, positivity=False, calc_llk=False, nb_cycle=1, show_cdi=False, fig_num=-1, zero_mask=False, confidence_interval_factor=0, confidence_interval_factor_mask_min=0.5, confidence_interval_factor_mask_max=1.2, update_psf=0, psf_filter=None, plot_axis=0)

RAAR cycle

Parameters:
  • positivity – apply a positivity restraint: the imaginary part inside the support is flipped.

  • calc_llk – if True, calculate llk while in Fourier space. If a positive integer is given, llk will be calculated every calc_llk cycle

  • nb_cycle – the number of cycles to perform

  • show_cdi – 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 intensity, as for ShowCDI()

  • zero_mask – if True, masked pixels (iobs<-1e19) are forced to zero, otherwise the calculated complex amplitude is kept with an optional scale factor.

  • confidence_interval_factor – a relaxation factor, with the projection of calculated amplitude being done towards the limit of the poisson confidence interval. A value of 1 corresponds to a 50% confidence interval, a value of 0 corresponds to a strict observed amplitude projection. [EXPERIMENTAL]

  • confidence_interval_factor_mask_max (confidence_interval_factor_mask_min,) – For masked pixels where a value has been estimated (e.g. with InterpIobsMask()), a confidence interval can be given as a factor to be applied to the interpolated observed intensity. This corresponds to values stored between -1e19 and -1e38. [EXPERIMENTAL]

  • update_psf – if >0, will update the partial coherence psf every update_psf cycles.

  • psf_filter – either None, “hann” or “tukey”: window type to filter the PSF update

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

op(cdi: CDI)

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).

class pynx.cdi.cl_operator.ScaleObj(method='I', verbose=False, lazy=False)

Scale the object according to the observed intensity. The scaling is either made against the amplitudes, the intensities, or the weighted intensities. This is only useful if a mask is used - the scale factor effectively only applies to masked intensities. :param method: ‘I’ (intensities), ‘F’ (amplitudes), ‘wI’ (weighted intensities), ‘P’ Poisson :return: nothing. The object is scaled to best match the intensities.

Parameters:
  • method – ‘I’ (intensities), ‘F’ (amplitudes), ‘wI’ (weighted intensities), ‘P’ (Poisson)

  • verbose – if True, print the scale factor

  • lazy – lazy evaluation

op(cdi)

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).

class pynx.cdi.cl_operator.ShowCDI(fig_num=None, i=None)
Parameters:
  • fig_num – the matplotlib figure number. if None, a new figure will be created each time.

  • i – if the object is 3D, display the ith plane (default: the center one)

static get_icalc(cdi, i=None)

This static, virtual function is used to get icalc, and should be derived depending on the GPU used.

Parameters:
  • cdi – the cdi object from which to extract the calculated intensity

  • i – if data is 3D, the index if the frame to extract

  • plot_axis – for 3D data, the plot_axis along which the cut plane will be selected

Returns:

the calculated intensity, as a float32 numpy array, or None if the intensity could not be calculated

class pynx.cdi.cl_operator.SupportExpand(n=1, update_nb_points_support=True)

Expand (or shrink) the support using a binary window convolution.

Parameters:
  • n – number of pixels to broaden the support, which will be done by a binary convolution with a window size equal to 2*n+1 along all dimensions. if n is negative, the support is instead shrunk, by performing the binary convolution and test on 1-support.

  • update_nb_points_support – if True (the default), the number of points in the support will be calculated and stored in the object

op(cdi)

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).

class pynx.cdi.cl_operator.SupportUpdate(threshold_relative=0.2, smooth_width=3, force_shrink=False, method='rms', post_expand=None, verbose=False, update_border_n=0, min_fraction=0, max_fraction=1, lazy=False)

Update the support

Update support.

Parameters:
  • threshold_relative – must be between 0 and 1. Only points with object amplitude above a value equal to relative_threshold * reference_value are kept in the support. reference_value can either: - use the fact that when converged, the square norm of the object is equal to the number of recorded photons (normalized Fourier Transform). Then: reference_value = sqrt((abs(obj)**2).sum()/nb_points_support) - or use threshold_percentile (see below, very slow, deprecated)

  • smooth_width – smooth the object amplitude using a gaussian of this width before calculating new support If this is a scalar, the smooth width is fixed to this value. If this is a 3-value tuple (or list or array), i.e. ‘smooth_width=2,0.5,600’, the smooth width will vary with the number of cycles recorded in the CDI object (as cdi.cycle), varying exponentially from the first to the second value over the number of cycles specified by the last value. With ‘smooth_width=a,b,nb’: - smooth_width = a * exp(-cdi.cycle/nb*log(b/a)) if cdi.cycle < nb - smooth_width = b if cdi.cycle >= nb

  • force_shrink – if True, the support can only shrink

  • method – either ‘max’ or ‘average’ or ‘rms’ (default), the threshold will be relative to either the maximum amplitude in the object, or the average or root-mean-square amplitude (computed inside support)

  • post_expand=1 – after the new support has been calculated, it can be processed using the SupportExpand operator, either one or multiple times, in order to ‘clean’ the support: - ‘post_expand=1’ will expand the support by 1 pixel - ‘post_expand=-1’ will shrink the support by 1 pixel - ‘post_expand=(-1,1)’ will shrink and then expand the support by 1 pixel - ‘post_expand=(-2,3)’ will shrink and then expand the support by respectively 2 and 3 pixels

  • verbose – if True, print number of points in support

  • update_border_n – if > 0, the only pixels affected by the support updated lie within +/- N pixels around the outer border of the support.

  • min_fraction – these are the minimum and maximum fraction of the support volume in the object. If the support volume fraction becomes smaller than min_fraction or larger than max_fraction, a corresponding exception will be raised. Example values: min_size=0.001, max_size=0.5

  • max_fraction – these are the minimum and maximum fraction of the support volume in the object. If the support volume fraction becomes smaller than min_fraction or larger than max_fraction, a corresponding exception will be raised. Example values: min_size=0.001, max_size=0.5

  • lazy – if True, this will be queued for later execution in the cdi object

Raises: SupportTooSmall or SupportTooLarge if support diverges according to {min|max}_fraction :returns: Nothing. self._support is updated

op(cdi)

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).

class pynx.cdi.cl_operator.UpdatePSF(nb_cycle=1, scale_in=1, filter=None)

Update the PSF while in detector space. Assumes that the calculated

Parameters:
  • nb_cycle – the number of cycle for the Richardson-Lucy update

  • scale_in – a scale factor by which the input values should be multiplied, typically because of FFT

  • filter – either None, “hann” or “tukey” - this will be used to filter the PSF update.

op(cdi: CDI)

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).

CDI Runner class

The ‘runner’ class is used for command-line analysis scripts.

class pynx.cdi.runner.runner.CDIRunner(argv, params, cdi_runner_scan_class)

Class to process CDI data, with parameters from the command-line or from a text file.

Parameters:
  • argv – the command-line parameters

  • params – parameters for the optimization, with some default values.

  • ptycho_runner_scan_class – the class to use to run the analysis.

check_mpi(task, timeout=60)
Parameters:
  • task – the string describing the task after which the check is made

  • timeout – the timeout in seconds after which the synchronisation will be considered failed, and an exception will be raised.

Returns:

nothing.

Raises:

a CDIRunnerException

check_params()

Check if self.params includes a minimal set of valid parameters If ‘help’ is given as command-line parameter, the help text is printed and the program exits.

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

Returns: Nothing. Will raise an exception if necessary

parse_arg(k, v=None)

Interprets one parameter. Will :param k: string with the name of the parameter :param v: value of the parameter, or None if a keyword parameter :return: True if parameter could be interpreted, False otherwise

parse_arg_beamline(k, v)

Parse argument in a beamline-specific way. This function only parses single arguments. If an argument is recognized and interpreted, the corresponding value is added to self.params

This method should be superseded in a beamline/instrument-specific child class.

Returns:

True if the argument is interpreted, false otherwise

parse_argv()

Parses all 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 as a python script :return: nothing. The parameters are accepted if understood, and stored in self.params

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.cdi.runner.runner.CDIRunnerException
class pynx.cdi.runner.runner.CDIRunnerScan(params, scan, timings=None)

Abstract class to handle CDI data. Must be derived to be used.

center_crop_data()

Auto-center diffraction pattern, and resize according to allowed FFT prime decomposition. Takes into account ROI if supplied. :return:

corr_flat_field()

Correct the image by a known flat field. Currently only applies correction for maxipix large pixels :return:

get_psf()

Get the PSF parameters from self.params[‘psf’]

Returns:

a tuple (model, fwhm, eta, update_psf) with the PSF parameters, or (psf, update_psf) if the PSF array was loaded from a file, or False if no psf is used.

init_mask()

Load mask if the corresponding parameter has been set, or just initialize an array of 0. This must be called after iobs has been imported, and before cropping, so at the end of load_data()

Returns:

Nothing. self.mask is created.

init_object()

Prepare the initial object. It will be random unless a file was supplied for input. :return: the initialised object, or an operator which can be used for a lazy initialisation

init_support()

Prepare the initial support. Note that the mask should be initialized first. if self.params[‘support’] == ‘auto’, nothing is done, and the support will be created after the cdi object, to use GPU functions. :return: Nothing. self.support is created

load_data()

Loads data. This function only loads data (2D or 3D) from generic files, and should be derived for beamline-specific imports (from spec+individual images, etc..)

Returns: nothing. Loads data in self.iobs, and initializes mask and initial support.

prepare_cdi()

Prepare CDI object from input parameters. If self.cdi already exists, it is re-used, and only the initial object and support are updated. To save memory, self.iobs and self.mask are set to None.

Returns:

nothing. Creates or updates self.cdi object.

prepare_data(init_mask=True, rebin_data=True, init_support=True, center_crop_data=True, corr_flat_field=True)

Prepare CDI data for processing.

Parameters:
  • init_mask – if True (the default), initialize the mask

  • rebin_data – if True (the default), rebin the data

  • init_support – if True (the default), initialize the support

  • center_crop_data – if True (the default), center & crop data

  • corr_flat_field – if True (the default), correct for the flat field

Returns:

prepare_processing_unit()

Prepare processing unit (CUDA, OpenCL, or CPU).

Returns: nothing. Creates self.processing_unit

rebin_data()

Will rebin the data (iobs and mask) and update pixel_size_detector according to self.params[‘rebin’]. this must be called at the end of load_data, after init_mask. If a ROI has been given, cropping will occur here. :return:

run(file_name=None, run_n=1)

Main work function. Will run selected algorithms according to parameters

Parameters:
  • file_name – output file_name. If None, the result is not saved.

  • run_n – the run number

Returns:

run_algorithm(algo_string, file_name=None, run_n=1)

Run a single or suite of algorithms in a given run.

Parameters:
  • algo_string – a single or suite of algorithm steps to use, which can either correspond to a change of parameters, i.e. ‘beta=0.5’, ‘support_threshold=0.3’, ‘positivity=1’, or operators which should be applied to the cdi object, such as: ‘hio**20’: 20 cycles of HIO ‘er’: a single cycle or error reduction ‘detwinhio**20’: 20 cycles of hio after halving the support ‘er**20*hio**100’: 100 cycles of hio followed by 20 cycles of ER ‘(sup*er**10*raar**100)**10’: 10 times [100 cycles of RAAR + 10 ER + Support update]

  • file_name – output file_name. If None, the result is not saved.

  • run_n – the run number

save_plot(file_name, algo_string=None)

Save a final plot for the object (only supported for 3D data) :param file_name: filename for the plot :param algo_string: will be used as subtitle :return: nothing

save_result(file_name=None, verbose=True)

Save the results (object, support,…) at the end of a run

Parameters:

file_name – the filename to save the data to. If its extension is either npz or cxi, this will override params[‘output_format’]. Otherwise, the extension will be replaced. Note that the full filename will be appended with the llk value.

Returns:

Regrid module reference

This module is used to transform far-field projections into a 3D interpolated dataset.

Rebuild the 3D reciprocal space by projecting a set of 2d speckle SAXS pattern taken at various rotation angles into a 3D regular volume

class pynx.cdi.runner.regrid.Regrid3D(mask, volume_shape, center, pixel_size, distance, slab_size=None, scale=None, ctx=None, devicetype='all', platformid=None, deviceid=None, block_size=None, memory=None, profile=False)

Project a 2D frame to a 3D volume taking into account the curvature of the Ewald’s sphere

Parameters:
  • mask – numpy array with the mask: needs to be of the same shape as the image

  • volume_shape – 3-tuple of int

  • center – 2-tuple of float (y,x)

  • pixel_size – float, size of the pixel in meter

  • distance – float, sample detector distance in meter

  • slab_size – Number of slices to be treated at one, the best is to leave the system guess

  • scale – zoom factor, 2 is like a 2x2x2 binning of the volume. Allows to work with smaller volumes.

  • ctx – actual working context, left to None for automatic initialization from device type or platformid/deviceid

  • devicetype – type of device, can be “CPU”, “GPU”, “ACC” or “ALL”

  • platformid – integer with the platform_identifier, as given by clinfo

  • deviceid – Integer with the device identifier, as given by clinfo

  • block_size – preferred workgroup size, may vary depending on the out come of the compilation

  • memory – minimum memory available on device

  • profile – switch on profiling to be able to profile at the kernel level, store profiling elements (makes code slightly slower)

calc_slabs()

Calculate the height of the slab depending on the device’s memory. The larger, the better

calc_slicing()

Calculate the slicing, i.e, for which slab in output, which lines of the image are needed

clean_slab()

Memset the slab

compile_kernels(kernel_files=None, compile_options=None)

Call the OpenCL compiler

Parameters:
  • kernel_files – list of path to the kernel (by default use the one declared in the class)

  • compile_options – string of compile options

get_slab()

After all frames have been projected onto the slab, retrieve it after normalization

Returns:

Ndarray of size (slab_size, volume_size_1, volume_size_2)

process_all(frames, oversampling_img=1, oversampling_rot=1)

Project all frames and rebuild the 3D volume

Parameters:
  • frames – dict with angle: frame as numpy.array

  • oversampling_img – Each pixel will be split in n x n and projected that many times

  • oversampling_rot – project multiple times each image between rot-d_rot/2 and rot+d_rot/2

Returns:

3D volume as numpy array

project_frames(l_frames, angles, step, vol_slice, img_slice=None, oversampling_img=1, oversampling_rot=1)

Project all frames onto the slab.

Parameters:
  • l_frames – list of frames

  • l_angles – angles associated with the frame

  • step – step size

  • vol_slice – fraction of the volume to use (slicing along Y!)

  • img_slice – fraction of the image to use

Returns:

the slab

project_one_frame(frame, rot, d_rot, vol_slice, img_slice=None, oversampling_img=1, oversampling_rot=1)

Projection of one image onto one slab :param frame: numpy.ndarray 2d, floa32 image :param rot: angle of rotation :param d_rot: angular step (used for oversampling_rot) :param vol_slice: Start/end row in the volume (slab along y) :param img_slice: Start/end row in the image :param oversampling_img: Each pixel will be split in n x n and projected that many times :param oversampling_rot: project multiple times each image between rot-d_rot/2 and rot+d_rot/2 :return: None

send_image(image, slice_=None)

Send image to the GPU

Parameters:
  • image – 2d numpy array

  • slice – slice object with the start and end of the buffer to be copied

Returns:

Nothing

send_mask(mask)

Send mask to the GPU

pynx.cdi.runner.regrid.parse_bliss_file(filename, title='dscan sz', rotation='ths', scan_len='1', callback=<function <lambda>>, maxi=None)

Scan a Bliss file and search for scans suitable for CXI image reconstruction

Parameters:
  • filname – str, name of the Bliss-Nexus file

  • title – the kind of scan one is looking for

  • rotation – name of the motor responsible for the rotation of the sample

  • scan_len – search only for scan of this length

  • callback – used for the progress-bar update

  • maxi – limit the search to this number of frames (used to speed-up reading in debug mode)

Returns:

dict with angle as key and image as value

Examples

# -*- coding: utf-8 -*-

# PyNX - Python tools for Nano-structures Crystallography
#   (c) 2017-present : ESRF-European Synchrotron Radiation Facility
#       authors:
#         Vincent Favre-Nicolin, favre@esrf.fr

import numpy as np
from scipy.fftpack import ifftshift, fftshift, fft2

# NB: It is possible to use the PYNX_PU environment variable to choose the GPU or language,
# e.g. using PYNX_PU=opencl or PYNX_PU=cuda or PYNX_PU=cpu or PYNX_PU=Titan, etc..

from pynx.utils.pattern import siemens_star
from pynx.cdi import *

# Note that this is a moderately difficult example: it is easy to get the outline of
# the siemens star, but the many voids within the structure can be difficult to
# reconstruct, especially when using the option to remove intensity in the centre.
# With this beamstop this needs to be run a number of time to get a really good result.

# Test on a simulated pattern (2D)
n = 512

# Siemens-Star object
obj0 = siemens_star(n, nb_rays=18, r_max=60, nb_rings=3, cheese_holes_nb=40)

# Start from a slightly loose disc support
x, y = np.meshgrid(np.arange(-n // 2, n // 2, dtype=np.float32), np.arange(-n // 2, n // 2, dtype=np.float32))
r = np.sqrt(x ** 2 + y ** 2)
support = r < 65

iobs = abs(ifftshift(fft2(fftshift(obj0.astype(np.complex64))))) ** 2
iobs = np.random.poisson(iobs * 1e10 / iobs.sum())
mask = np.zeros_like(iobs, dtype=np.int16)
if True:
    # Mask some values in the central beam
    print("Removing %6.3f%% intensity (more difficult, retry to get a good result" %
          (iobs[255:257, 255:257].sum() / iobs.sum() * 100))
    iobs[255:257, 255:257] = 0
    mask[255:257, 255:257] = 1

cdi = CDI(fftshift(iobs), obj=None, support=fftshift(support), mask=fftshift(mask), wavelength=1e-10,
          pixel_size_detector=55e-10)

cdi.init_free_pixels()

# Init real object from the chosen support
cdi = InitObjRandom(src="support", amin=0.8, amax=1, phirange=0) * cdi

# Initial scaling of the object [ only useful if there are masked pixels !]
cdi = ScaleObj(method='F') * cdi

show = 40

# Do 100 cycles of RAAR
cdi = RAAR(calc_llk=20) ** 100 * cdi

if show > 0:
    cdi = ShowCDI(fig_num=1) * cdi

sup = SupportUpdate(threshold_relative=0.3, smooth_width=(2.0, 0.5, 800), force_shrink=False)

cdi = (sup * ER(calc_llk=20, show_cdi=show, fig_num=1) ** 5 * RAAR(calc_llk=20, show_cdi=show,
                                                                   fig_num=1) ** 40) ** 5 * cdi
cdi = DetwinRAAR(nb_cycle=10) * cdi
cdi = (sup * ER(calc_llk=20, show_cdi=show, fig_num=1) ** 5 * RAAR(calc_llk=20, show_cdi=show,
                                                                   fig_num=1) ** 40) ** 20 * cdi

# Finish with ML or ER
# cdi = ML(reg_fac=1e-2, calc_llk=5, show_cdi=show, fig_num=1) ** 100 * cdi
cdi = ER(calc_llk=20, show_cdi=show, fig_num=1) ** 100 * cdi

cdi = ShowCDI(fig_num=1) * cdi
# -*- coding: utf-8 -*-

# PyNX - Python tools for Nano-structures Crystallography
#   (c) 2017-present : ESRF-European Synchrotron Radiation Facility
#       authors:
#         Vincent Favre-Nicolin, favre@esrf.fr

import numpy as np
from scipy.ndimage import gaussian_filter
from scipy.fftpack import fftshift
import fabio
from matplotlib import pyplot as plt

from pynx.cdi import *

# CDI example from an experimental data set from id10 (courtesy of Yuriy Chushkin)

iobs = np.flipud(fabio.open("data/logo5mu_20sec.edf").data)
support = np.flipud(fabio.open("data/mask_logo5mu_20sec.edf").data)
n = len(iobs)
x, y = np.meshgrid(np.arange(0, n, dtype=np.float32), np.arange(0, n, dtype=np.float32))

# Mask specific to this dataset (from beamstop, after symmetrization of observed data)
mask = np.logical_or(iobs != 0, np.logical_or(abs(x - 258) > 30, abs(y - 258) > 30))
mask *= np.logical_or(iobs != 0, np.logical_or(abs(x - 246) > 30, abs(y) < 495))
mask *= np.logical_or(iobs != 0, np.logical_or(abs(x - 266) > 30, abs(y) > 20))
mask *= np.logical_or(iobs != 0, np.logical_or(abs(x - 10) > 30, abs(y - 270) > 5))
mask *= np.logical_or(iobs != 0, np.logical_or(abs(x - 498) > 30, abs(y - 240) > 5))
mask *= np.logical_or(iobs != 0, np.logical_or(abs(x) > 30, abs(y) > 20))
mask *= np.logical_or(iobs != 0, np.logical_or(abs(x - 510) > 30, abs(y - 510) > 20))
mask = (mask == 0)  # 0: OK, 1: masked

plt.figure(1, figsize=(8, 8))

# ========================= Try first from the known support (easy !) =============================
cdi = CDI(fftshift(iobs), obj=None, support=fftshift(support), mask=fftshift(mask), wavelength=1e-10,
          pixel_size_detector=55e-6)

# Set free mask for unbiased likelihood estimation
cdi.init_free_pixels()

# Init real object from the supplied support
cdi = InitObjRandom(src="support", amin=0.8, amax=1, phirange=0) * cdi

# Initial scaling, required by mask
cdi = ScaleObj(method='F') * cdi

# Do 4 * (50 cycles of HIO + 20 of ER), displaying object after each group of cycle
cdi = (ShowCDI(fig_num=1) * ER(calc_llk=10) ** 20 * HIO(calc_llk=20) ** 50) ** 4 * cdi

# cdi = (ShowCDI(fig_num=1) * ER(calc_llk=10) ** 20 * CF(calc_llk=20) ** 50) ** 4 * cdi

# cdi = (ShowCDI(fig_num=1) * ER(calc_llk=10) ** 20 * RAAR(calc_llk=20) ** 50) ** 4 * cdi

print("\n======================================================================\n")
print("This was too easy - start again from a loose support !\n")

# ========================= Try again from a loose support ========================================
support = np.flipud(fabio.open("data/mask_logo5mu_20sec.edf").data)
support = gaussian_filter(support.astype(np.float32), 4) > 0.2
cdi = CDI(fftshift(iobs), obj=None, support=fftshift(support), mask=fftshift(mask), wavelength=1e-10,
          pixel_size_detector=55e-6)

# Set free mask for unbiased likelihood estimation
cdi.init_free_pixels()

# Init real object from the chosen support
cdi = InitObjRandom(src="support", amin=0.8, amax=1, phirange=0) * cdi

# Initial scaling, required by mask
cdi = ScaleObj(method='F') * cdi

# Display every N cycles
show = 100

# Do 200 cycles of HIO, displaying object every N cycle and log-likelihood every 20 cycle
cdi = HIO(calc_llk=20, show_cdi=show, fig_num=2) ** 200 * cdi

# Support update operator
sup = SupportUpdate(threshold_relative=0.25, smooth_width=(2, 0.5, 800), force_shrink=False)

if True:
    # Do 40 cycles of HIO, then 5 of ER, update support, repeat
    cdi = (sup * ER(calc_llk=20, show_cdi=show, fig_num=2) ** 5
           * HIO(calc_llk=20, show_cdi=show, fig_num=2) ** 40) ** 20 * cdi
else:
    # Do 40 cycles of HIO, update support, repeat
    cdi = (sup * HIO(calc_llk=20, show_cdi=show, fig_num=2) ** 40) ** 20 * cdi

# Finish with ML or ER
# cdi = ML(reg_fac=1e-2, calc_llk=20, show_cdi=show, fig_num=1) ** 100 * cdi
cdi = ER(calc_llk=20, show_cdi=show, fig_num=2) ** 100 * cdi