PyDDM.ddm_analysis_and_fitting.DDM_Analysis#
- class DDM_Analysis(data_yaml, load_images=True)#
Bases:
object
Performs preprossing of data, such as cropping, windowing etc. DDM calculations are performed on the processed time series to produce a xarray DataSet with DDM matrix, radial averages and ISF. The analysis parameters are provided by the user in a YAML file: More information here
Methods
Calculates the DDM matrix This function computes the DDM matrix.
Create two-time correlation matrix
Generates plot, which can be saved as PDF, to summarize the DDM matrix calulations.
Opens yaml file and extracts parameter values.
For more info, see Colin, R., Zhang, R.
Resave the DDM dataset.
Change file name to save the data to disk.
Based off user-provided parameters, prepares images for the DDM analysis.
Creates the xarray Dataset and PDF report.
- calculate_DDM_matrix(quiet=False, velocity=[0, 0], bg_subtract_for_AB_determination=None, **kwargs)#
Calculates the DDM matrix This function computes the DDM matrix. The radially averaged DDM matrix will then also be found, along with estimates of the background and amplitude. From the amplitude and background, we can extract the intermediate scattering function (ISF) from the DDM matrix. All of these computed variables will be stored as an xarray Dataset.
- Parameters:
quiet (boolean (optional)) – If set to False, then when calculating the DDM matrix, a message will print out at about every fourth time lag calculated (with a timestamp)
velocity (array-like (optional)) – Velocity in x and y direction.
bg_subtract_for_AB_determination ({None, 'mode', 'median'}) – Deafault option is None. In method to find A(q) and B, one can subtract the mode or median of the images.
**overlap_method ({0,1,2,3}, optional) – Optional keyword argument. Will be set to 2 if not specified here nor in the YAML file. Determines how overlapped the different image pairs are. Let’s say you are finding all pairs of images separated by a lag time of 10 frames. You’d have frame 1 and 11, 2 and 12, 3 and 13, etc. One could use each possible pair when calculating the DDM matrix (overlap_method=3, maximally overlapping pairs). One could only look at non-overlapping pairs, like frame 1 and 11, 11 and 21, 21 and 31, etc (overlap_method=0, non-overlapping). Or one could do something in between. For overlap_method=2 (the DEFAULT), there will be about 3 or 4 pairs (at most) between two frames. That is, if we have a lag time of 10 frames, we will use the pairs of frame 1 and 11, 5 and 15, 9 and 19, etc. For overlap_method=1, we do something similar to overlap_method=2, but we compute at most number_differences_max (default:300) image differences per lag time. So for example, with a lag time of 1 frame and a movie that has 1000 frames, we could theoretically use 999 differences of images (and we would for the other methods). But for overlap_method=2, we would only use 50 of those. This is for quickening up the computation.
**background_method ({0,1,2,3,4,5}, optional) – Optional keyword argument. Will be set to 0 if not specified here nor in the YAML file. Determines how we estimate the \(B\) parameter. This can be done by lookinag the average of the Fourier transform of each image (not image difference!) squared (that’s background_method=0). We could also use the minimum of the DDM matrix (background_method=1). Or we could use the average DDM matrix for the largest q (background_method=2). Or we could just set \(B\) to zero (background_method=3).
**amplitude_method ({0,1,2}, optional) – Optional keyword argument.
**number_lag_times (int) – Optional keyword argument. Must be set in the YAML file. You may pass this optional keyword argument if you want to overwrite the value for the number of lag times set in the YAML file.
- Returns:
ddm_dataset – Dataset containing the DDM matrix, the estimate of the background and amplitude, the ISF, etc. Coordinates of these variables will be the wavevector (either as q_x, q_y or the magnitude q), lagtime, etc.
- Return type:
xarray Dataset
- createTwoTimeCorr(ddm_var_dataset, qindex)#
Create two-time correlation matrix
After generating the DDM matrix as a function of time lag as well as of time, can use this function to generate a 2D two-time correlation function for a particular wavenumber.
- Parameters:
ddm_variability (xarray dataset) – Results of function ‘variationInDDMMatrix’ in ‘ddm_analysis_and_fitting’ code
q_index (int) – Index of the wavenumber array
- Returns:
twotimecorrelation – Two time correlation function
- Return type:
array
- find_alignment_factor(ddmmatrix3d, orientation_axis=0, remove_vert_line=True, remove_hor_line=True)#
- Parameters:
orientation_axis (TYPE, optional) – DESCRIPTION. The default is np.pi/4.
- Return type:
None.
- find_alignment_factor_one_lagtime(ddmmatrix2d, orientation_axis=0, remove_vert_line=True, remove_hor_line=True)#
- Parameters:
orientation_axis (TYPE, optional) – DESCRIPTION. The default is np.pi/4.
- Return type:
None.
- generate_plots(ddmdataset, pdf_to_save_to=None, q_to_see=1.5, num=None)#
Generates plot, which can be saved as PDF, to summarize the DDM matrix calulations. The region of interest is displayed; a Fourier transform of a difference image; a graph of radial averages of FFT of frames, used to determine amplitude and background and a plot of the intermediate scattering function at one given q-value.
- Parameters:
ddmdataset (Dataset) – DDM Dataset
pdf_to_save_to (TYPE, optional)
q_to_see (float) – The magnitude of the wavector for which to plot the ISF units: μm$^{-1}$
num (int, optional) – The number of the ROI in case multiple ROIs have been analyzed
- q_to_see: float, optional
q value to look at ISF
- loadYAML()#
Opens yaml file and extracts parameter values.
To analyze a recorded movie, there must be associated metadata. That can either be in the form of a yaml file or string in the yaml format. The yaml file
- phiDM(lagt, halfsize, use_gf=True, gfsize=3, err_limit=2e-05)#
For more info, see Colin, R., Zhang, R. & Wilson, L. G. Fast, high-throughput measurement of collective behaviour in a bacterial population. Journal of The Royal Society Interface 11, 20140486 (2014). https://royalsocietypublishing.org/doi/10.1098/rsif.2014.0486
- Parameters:
lagt (int) – Lag time (in units of frames).
halfsize (int) – Centered above zero frequency, size of region to fit to plane.
use_gf (boolean, optional) – Use Gaussian filter or not. The default is True.
gfsize (int, optional) – Size of Gaussian filter to apply. The default is 3.
err_limit (float, optional) – If error between data and plane that is fit is greater than this value per pixel, then won’t be used to get avg velocity
- Returns:
phiDM_dataset – Dataset contains the phase of the Fourier tranform of each image, the velocity (v_x and v_y), and error between the fit of plane and the difference in phase.
- Return type:
xarray Dataset
- resave_ddm_dataset(ddmdataset, file_name=None)#
Resave the DDM dataset.
- Parameters:
ddmdataset (xarray dataset) – DDM dataset
file_name (string or None, optional) – Filename for saving data. The default is None. If None, then will use what is set in the attribute ‘filename_for_saving_data’. Note that ‘_ddmmatrix.nc’ is appended to the end of the filename.
- set_filename_for_saving(filename, quiet=False)#
Change file name to save the data to disk. This is the file that will store the ddm_dataset as a netCDF file (with extension .nc).
- Parameters:
filename (str) – New file name
- setup(load_images)#
Based off user-provided parameters, prepares images for the DDM analysis.
Using the parameters under the ‘Analysis_parameters’ section of the YAML file, this will get the image stack ready for DDM analysis. Possible things to do done include cropping the image to focus on a particular region of interest (ROI), binnning the image, applying a windowing function.
- variationInDDMMatrix(lagtime, orientation_axis=0, save_full_ddmmat=True, velocity=[0, 0])#
Creates the xarray Dataset and PDF report.
- Parameters:
lagtime (int or list-like) – If an integer, will compute the DDM matrix for this lag time over all times. If an arrary or list, then all lag times in that array/list will be used.
orientation_axis (float, optional) – Axis for computing the alignment factor of the 2D DDM matrix for a given lag time and time
save_full_ddmmat (bool, optional) – If True, then will save the DDM matrix as a function of q_x and q_y for each lag time and time. If False, then it will radially average that matrix so that the DDM matrix is just a function of the magnitude of the wavevector, the lagtime, and the time. If True, then it may potentially use up a lot of memory.
velocity (list-like, optional) – Deafult is [0,0]. If not [0,0], then will use phiDM to correct for drift or ballistic motion.
- Returns:
ddm_dataset – Dataset containing the DDM matrix and associated data and metadata.
- Return type:
xarray Dataset