PointSpreadFunction calculation

psf(::Type{ModeWidefield}, sz::NTuple, pp::PSFParams; sampling=get_sampling(sz, pp))

Calculates the widefield single-frequency point spread function (psf), i.e. the image of a single (very small) emitter. Most of the parameters (such as refractive index, numerical aperture, vacuum wavelength, aberrations etc.) are hidden in the parameter structure argument pp, which should be generated via the PSFParams() constructor. See `PSFParams() for details.

Parameters

• sz: size tuple of the final PSF
• pp: PSF parameters of the PSF. See PSFParams() for details. The argument pp.aplanatic defined whether an excitation or emission PSF is calculated.
• sampling=nothing: The sampling parameters of the resulting PSF. If nothing is provided, the PSF will be sampled according to the Abbe limit.
• use_resampling=true: Exploits a calculation trick, which first calculates the amplitudes are calculated on a twice coarser grid and the result is upsampled. This increases the speed but may be less accurate.
• return_amp=false: Has to be false since confocal amplitude spread functions do not exist for non-zero pinhole sizes.

See also:

• apsf: calculates the underlying amplitude point spread function (apsf)

Example:

julia> pp = PSFParams(0.5,1.4,1.52); sz=(128,128,128); sampling=(0.050,0.050,0.200);
# an emission PSF of an isotropic (freely rotating emitter)

julia> p_wf = psf(sz, pp; sampling=sampling);
# an emission PSF of an emitter oriented along the Z direction

julia> pp_dipole = PSFParams(0.5,1.4,1.52; transition_dipole=[0.0,0.0,1.0]);

julia> p_dipole = psf(sz, pp_dipole; sampling=sampling);

psf(::Type{ModeConfocal}, sz::NTuple, pp_em::PSFParams; pp_ex=nothing, pinhole=nothing, pinhole_ft=disc_pinhole_ft, sampling=nothing, use_resampling=true, return_amp=false, pinhole_positions=[(0.0,0.0)], ex_modifier=modify_ident) # unclear why the resampling seems to be so bad

Calculates a confocal point spread function. The normalisation is such that a completely open pinhole diameter yields the excitation PSF with its normalization. Returns the PSF or a vector of PSFs.

Parameters

• sz: size tuple of the final PSF
• pp_em: PSF parameters of the emission PSF. PSF parameters of the PSF. See PSFParams() for details. This should include the emission wavelength
• pp_ex=nothing: This is a required named parameter, containing all the settings for the excitation PSF. This should include the excitation wavelength as well as typically aplanatic=aplanatic_illumination.
• pinhole=nothing: The diameter of the pinhole in Airy Units (AU = 1.22 λ/NA). A pinhole size of one AU corresponds to a pinhole border falling onto the first zero of a corresponding paraxial emission PSF.
• sampling=nothing: The sampling parameters of the resulting PSF.
• use_resampling=true: Exploits a calculation trick, which first calculates the individual PSFs on a twice coarser grid in XY and Z and then upsamples the result. Note that this may be inappropriate due to undersampling due to the Stokes shift which is neglected here. But warnings will then result during the calculations of the subsampled widefield PSFs.
• return_amp=false: Has to be false since confocal amplitude spread functions do not exist for non-zero pinhole sizes.
• pinhole_positions=[(0.0,0.0)]: A list of pinhole positions in pixels. One PSF will be returned for each pinhole position. If only a single pinhole position is supplied the PSF will directly be returned instead of a vector of PSFs. Specifies the precise position(s) of the pinholes in the detection path. This allows to simulate an offset (misadjusted) pinhole, or (as a vector of tuples) a PSF for a whole set of positions. See the psf(ModeISM, ...) for more details.
• pinhole_ft=disc_pinhole_ft: Specifies which function is used to calculate the Fourier transform of the pinhole. This allows the user to control the pinhole shape.
• ex_modifier: A function that can modify the excitation PSF. By default the identity is used, but other options are modify_square to calculate a two-photon confocal PSF. However, you should typically use the Mode2Photon to do this. Note that the order of the PSFParams are reversed. By default (modify_ident), no modication is performed.

julia> pp_em = PSFParams(0.5,1.4,1.52; mode=ModeConfocal);
julia> pp_ex = PSFParams(pp_em; λ=0.488, aplanatic=aplanatic_illumination);
julia> p_conf = psf((128,128,128),pp_em; pp_ex=pp_ex, pinhole=0.1, sampling=(0.040,0.040,0.100));

psf(::Type{ModeLightsheet}, sz::NTuple, pp_em::PSFParams; pp_ex=nothing, sampling=nothing, use_resampling=true, dynamically_scanned=false) # unclear why the resampling seems to be so bad

Calculates a point spread function for a lightsheet microscopy. Not that the coordinate system refers to the detection PSF. The NA of the excitation PSF refers to the Z and Y direction. The polarisation of the excitation PSF along X means along Z and along Y means along Y in the final volume. Note that the excitation is assumed for the best focus position (beam waist). In reality the PSF will be spatially dependent due to the excitation focus.

Parameters

• sz: size tuple of the final PSF
• pp_em: PSF parameters of the emission PSF. PSF parameters of the PSF. See PSFParams() for details. This should include the emission wavelength
• pp_ex=nothing: This is a required named parameter, containing all the settings for the excitation PSF. This should include the excitation wavelength as well as typically aplanatic=aplanatic_illumination.
• pinhole=nothing: The diameter of the pinhole in Airy Units (AU = 1.22 λ/NA). A pinhole size of one AU corresponds to a pinhole border falling onto the first zero of a corresponding paraxial emission PSF.
• sampling=nothing: The sampling parameters of the resulting PSF.
• use_resampling=true: Exploits a calculation trick, which first calculates the individual PSFs on a twice coarser grid in XY and Z and then upsamples the result. Note that this may be inappropriate due to undersampling due to the Stokes shift which is neglected here. But warnings will then result during the calculations of the subsampled widefield PSFs.
• 'dynamically_scanned': specifies, whether the lightsheet uses a dynamically scanned beam. This is important for the calculation of the PSF. If the beam is scanned, the PSF will be calculated in the same way as for a confocal PSF. If the beam is not scanned, the excitation PSF will be assuming a focussing via a cylinder lens.
• ex_modifier: A function that can modify the excitation PSF. By default the identity is used, but other options are modify_square to calculate a two-photon confocal PSF. However, you should typically use the Mode2Photon to do this. Note that the order of the PSFParams are reversed. By default (modify_ident), no modication is performed.
• sigmaz = nothing: If the user provides a sigmaz value, the excitation calculation is using a Guassian excitation PSF with the specified sigma. pp_ex is ignored.

julia> pp_em = PSFParams(0.5,0.9,1.33; mode=ModeLightsheet);
julia> pp_ex = PSFParams(0.488, 0.2, 1.33);
julia> p_lightsheet = psf((128,128,128), pp_em; pp_ex=pp_ex, sampling=(0.040,0.040,0.100), dynamically_scanned=true);
julia> p_lightsheet_gauss = psf((128,128,128), pp_em; sampling=(0.040,0.040,0.100), sigma_z=1.0); # Gaussian excitation PSF

psf(::Type{ModeISM}, sz::NTuple, pp_em::PSFParams; pinhole=nothing, sampling=nothing, pinhole_ft=box_pinhole_ft, pinhole_positions=[(0.0,0.0)], pinhole_dist=0.5, pinhole_grid=(5,5), ism_pos=ism_positions_rect, args...)

Calculates a confocal point spread function. The normalisation is such that a completely open pinhole diameter yields the excitation PSF with its normalization. Returns the PSF or a vector of PSFs.

Parameters

• sz: size tuple of the final PSF
• pp_em: PSF parameters of the emission PSF. This should include the emission wavelength
• pp_ex=nothing: This is a required named parameter, containing all the settings for the excitation PSF. This should include the excitation wavelength as well as typically aplanatic=aplanatic_illumination.
• pinhole=nothing: The diameter of each pinhole in Airy Units (AU = 1.22 λ/NA). By default (pinholes=nothing) the pinhole size is automatically calculated by the pinhole_dist parameter below to yield mutually touching pinholes.
• sampling=nothing: The sampling parameters of the resulting PSF.
• use_resampling=true: Exploits a calculation trick, which first calculates the individual PSFs on a twice coarser grid in XY and Z and then upsamples the result. Note that this may be inappropriate due toundersampling due to the Stokes shift which is neglected here. But warnings will then result during the calculations of the subsampled widefield PSFs.
• return_amp=false: Has to be false since confocal amplitude spread functions do not exist for non-zero pinhole sizes.
• pinhole_positions=nothing: A list of pinhole positions in pixel coordinates. One PSF will be returned for each pinhole position. If only a single pinhole position is supplied the PSF will directly be returned instead of a vector of PSFs. Be careful, since this is not the same as the position of the relative shift of the images.
• pinhole_dist=0.5: a value or tuple specified the distances between pinholes, when arranged in a grid.
• pinhole_ft=box_pinhole_ft: Specifies which function is used to calculate the Fourier transform of the pinhole. This allows the user to control the pinhole shape. E.g. hexagonal pattern with round pinholes

julia> sz=(128,128,128); sampling = (0.04,0.04,0.120)
julia> pp_em = PSFParams(0.5,1.4,1.52; mode=ModeISM);
julia> pp_ex = PSFParams(pp_em; λ=0.488, aplanatic=aplanatic_illumination);
julia> p_ism = psf(sz,pp_em; pp_ex=pp_ex, pinhole=0.21, pinhole_dist=0.2, sampling=sampling);

psf(::Type{Mode2Photon}, sz::NTuple, pp_ex::PSFParams; pp_em=nothing, pinhole=nothing, sampling=nothing, pinhole_ft=disc_pinhole_ft, args...)

Calculates a 2-photon (potentially confocal) point spread function. Returns the PSF or a vector of PSFs.

Parameters

• sz: size tuple of the final PSF
• pp_ex: PSF parameters of the excitation PSF. This should include the exission wavelength (typically in the IR region). Please make sure to also set pp.aplanatic=aplanatic_illumination.
• pp_em: PSF parameters of the emission PSF. This only needs to be supplied, if a pinhole is used. Otherwise non-descanned detection (NDD) is assumed.
• pinhole=nothing: If nothing, NDD is assumed and the PSF is only the square of the excitation PSF. The diameter of each pinhole in Airy Units (AU = 1.22 λ/NA).
• sampling=nothing: The sampling parameters of the resulting PSF.
• use_resampling=true: Exploits a calculation trick, which first calculates the individual PSFs on a twice coarser grid in XY and Z and then upsamples the result. Note that this may be inappropriate due toundersampling due to the Stokes shift which is neglected here. But warnings will then result during the calculations of the subsampled widefield PSFs.
• return_amp=false: Has to be false since confocal amplitude spread functions do not exist for non-zero pinhole sizes.
• pinhole_positions=[(0.0,0.0)]: A list of pinhole positions in pixels. One PSF will be returned for each pinhole position. If only a single pinhole position is supplied the PSF will directly be returned instead of a vector of PSFs.

julia> sz=(128,128,128); sampling = (0.04,0.04,0.120)
julia> pp_ex = PSFParams(0.8,1.4,1.52; mode=Mode2Photon, aplanatic= aplanatic_illumination);
julia> p_2p = psf(sz,pp_ex; sampling=sampling);

psf(::Type{ModeSTED}, sz::NTuple, pp_ex::PSFParams, pp_sted::PSFParams; pp_em=nothing, pinhole=nothing, sampling=nothing, pinhole_ft=disc_pinhole_ft, args...)

Calculates a 2-photon (potentially confocal) point spread function. Returns the PSF or a vector of PSFs.

Parameters

• sz: size tuple of the final PSF
• pp_ex: PSF parameters of the excitation PSF. This should include the exission wavelength (typically in the IR region). Please make sure to also set pp.aplanatic=aplanatic_illumination.
• pp_sted: PSF parameters of the excitation PSF. This should include the exission wavelength (typically in the IR region). Please make sure to also set pp.aplanatic=aplanatic_illumination.
• pp_em: PSF parameters of the emission PSF. This only needs to be supplied, if a pinhole is used. Otherwise non-descanned detection (NDD) is assumed.
• pinhole=nothing: If nothing, NDD is assumed and the PSF is only the square of the excitation PSF. The diameter of each pinhole in Airy Units (AU = 1.22 λ/NA).
• sampling=nothing: The sampling parameters of the resulting PSF.
• use_resampling=true: Exploits a calculation trick, which first calculates the individual PSFs on a twice coarser grid in XY and Z and then upsamples the result. Note that this may be inappropriate due toundersampling due to the Stokes shift which is neglected here. But warnings will then result during the calculations of the subsampled widefield PSFs.
• return_amp=false: Has to be false since confocal amplitude spread functions do not exist for non-zero pinhole sizes.
• pinhole_positions=[(0.0,0.0)]: A list of pinhole positions in pixels. One PSF will be returned for each pinhole position. If only a single pinhole position is supplied the PSF will directly be returned instead of a vector of PSFs.

julia> sz=(128,128,128); sampling = (0.04,0.04,0.120)
julia> pp_em = PSFParams(0.510,1.4,1.52; mode=ModeSTED); # STED suppression beam
julia> pp_ex = PSFParams(0.488,1.4,1.52; aplanatic= aplanatic_illumination);
julia> pp_sted = PSFParams(0.570,1.4,1.52; method=MethodPropagateIterative, aplanatic= aplanatic_illumination, pol=pol_circ_spiral);
julia> p_sted_01 = psf(sz, pp_em; pp_ex=pp_ex, pinhole= 2.0, sampling=sampling, pp_sted=pp_sted, rel_sted_intensity=0.1);
julia> p_sted_1 = psf(sz, pp_em; pp_ex=pp_ex, pinhole= 2.0, sampling=sampling, pp_sted=pp_sted, rel_sted_intensity=1.0);
julia> p_sted_5 = psf(sz, pp_em; pp_ex=pp_ex, pinhole= 2.0, sampling=sampling, pp_sted=pp_sted, rel_sted_intensity=5.0);
julia> pp_sted_b = PSFParams(0.570,1.4,1.52; method=MethodPropagateIterative, aplanatic= aplanatic_illumination, pol=pol_circ_tophat);
julia> p_sted_b5 = psf(sz, pp_em; pp_ex=pp_ex, pinhole= 2.0, sampling=sampling, pp_sted=pp_sted_b, rel_sted_intensity=5.0);
julia> using View5D; @vt p_sted_01 p_sted_1 p_sted_5 p_sted_b5

psf(::Type{Mode4Pi}, sz::NTuple, pp_ex::PSFParams; sampling=nothing, pp_ex2=pp_ex, pp_em=nothing, pp_em2=nothing, rel_ex_phase = 0.0, rel_em_phase = 0.0, pinhole=nothing, pinhole_ft=disc_pinhole_ft, pinhole_positions=[(0.0,0.0)], ex_modifier=modify_square)

Calculates a 4Pi point spread function. Note that the default is using a two-photon excitation. Returns the PSF or a vector of PSFs.

Parameters

• sz: size tuple of the final PSF
• pp_ex: PSF parameters of the excitation PSF. This should include the exission wavelength (typically in the IR region). Please make sure to also set pp.aplanatic=aplanatic_illumination.
• pp_ex2: PSF parameters of the other side excitation PSF. If nothing is provided, single-sided excitation (e.g. 4Pi Type B) is assumed.
• ex_modifier=modify_square: If modify_square, two-photon excitation is assumed. Use modify_ident for single-photon excitation.
• pp_em: PSF parameters of the emission PSF. This only needs to be supplied, if a pinhole is used.
• pp_em2: PSF parameters of the other side emission PSF. If nothing is supplied, Type A 4Pi microscopy is assumed with single-sided (or incoherent) detection.
• pinhole=nothing: If nothing, NDD is assumed and the PSF is only the square of the excitation PSF. The diameter of each pinhole in Airy Units (AU = 1.22 λ/NA).
• sampling=nothing: The sampling parameters of the resulting PSF.
• pinhole_positions=[(0.0,0.0)]: A list of pinhole positions in pixels. One PSF will be returned for each pinhole position. If only a single pinhole position is supplied the PSF will directly be returned instead of a vector of PSFs.

julia> sampling = (0.04,0.04,0.04)
julia> sz = (128,128,128)
julia> pp_em = PSFParams(0.5,1.3,1.52; mode=Mode4Pi, pol=pol_x);
julia> pp_ex = PSFParams(0.800,1.3,1.52; aplanatic=aplanatic_illumination, mode=Mode4Pi, pol=pol_x);
julia> @time p_4pi = psf(sz,pp_ex; pp_ex2=pp_ex, pp_em=pp_em, pp_em2=pp_em, sampling=sampling, pinhole=2.0, ex_modifier=modify_square);

psf(sz::NTuple, pp::PSFParams; sampling=get_sampling(sz, pp))

Calculates the point spread function (psf), i.e. the image of a single (very small) emitter. Most of the parameters (such as refractive index, numerical aperture, vacuum wavelength, aberrations etc.) are hidden in the parameter structure argument pp, which should be generated via the PSFParams() constructor. See `PSFParams() for details. Note that the field pp.mode defines the microscopic mode to simulate. Currently implemented are the default ModeWidefield and ModeConfocal.

See also:

• apsf(): calculates the underlying amplitude point spread function (apsf)

Example:

julia> pp = PSFParams(0.5,1.4,1.52);
julia> p = psf((128,128,128),pp; sampling=(0.050,0.050,0.200));

apsf(::Type{MethodParaxial}, sz::NTuple, pp::PSFParams; sampling=nothing) 
Calculates a paraxial amplitude PSF. Typically `pp.polarization` should be `pol_scalar`. However, other polarisation types yield two channels in the 4th dimension.
One for X and one for Y polarisation. In the paraxial approximation there is no Z polarisation.

apsf(::Type{MethodCZT}, sz::NTuple, pp::PSFParams; sampling=nothing, center_kz=false) 
Calculate amplitude PSF using Chrip-Z transform. the PSF at a given z_depth requires a bigger window to avoid or reduce wrap-around effect of FFT operation.
Firstly given an input parameter z_depth to calculate the xy_plane zoom factor c then apply on the pupil based on the desired depth. Calculate the inverse CZT of the 
field at pupil and zoom-out by the factor c.
+ z_depth: Z-Position from the focus to propagate to. If calcualte 3D apsf (sz[3] is larger than one), a stack is generated with z_depth being the position of the middle (floor.(size/2)) slice.
+ c_want: plane zoom factor only vary by λ, NA, sampling and z_depth. c_allow = λ/(NA*sampling)

apsf(sz::NTuple, pp::PSFParams; sampling=nothing, center_kz=false)

dispatches to various amplitude point spread function calculation routines. Note that the method entry in pp defines which calculation method to be used. Alternatively a different method can be chosen like this: apsf(::Type{MethodShell}, sz, ..).

Arguments:

• sz: NTuple of size to generate
• pp: PSF parameter structure, also contains the dtype. See PSFParam() for details
• sampling: NTuple for pixel pitch information
• center_kz: if true, the McCutchen pupil will be centered along the kz direction. This is important to be able to apply a consecutive resampling without errors. However, the phase values are then not correct, which does not matter for intensity PSFs.

See also:

• psf(): calculates the intensity point spread function (psf) by taking (sum along field components of) the absolute square of the corresponding apsf.

Example:

julia> pp = PSFParams(500.0,1.4,1.52)
julia> p = apsf((128,128,128),pp; sampling=(50,50,100));

modify_ident(h)

A function from a series of modficator functions which serve to modify the excitation intensity. This particular version does nothing. It just returns the input h.

modify_square(h)

A function from a series of modficator functions which serve to modify the excitation intensity. This version returns the absolute square of the input h. It is useful for two-photon microscopy.

modify_STED_exp(h, h_doughnut)

A function from a series of modficator functions which serve to modify the excitation intensity. This version applies a second STED PSF h_doughnut to deplete the excitation PSF h. The depletion is modeled by an exponential decay as obtained by a relatively short STED pulse following a very short excitation pulse. Note the h_doughnut should be the depletion PSF. Note that this function should be used as a tool to define a new function with fixed second argument.

modify_STED_hyper(h, h_doughnut)

A function from a series of modficator functions which serve to modify the excitation intensity. This version applies a second STED psf h_doughnut to deplete the excitation psf h. The depletion is modeled by an hyperbolic saturation curve based on anaylsing the steady-state kinetics. This is mostly suitable for CW-STED. Note the h_doughnut should be the depletion PSF. Note that this function should be used as a tool to define a new function with fixed second argument. The normalizations are such that h as well as h_doughnut reache 50% saturation at h==1.0. This means that typically your want to use max(h) << 1 and max(h_doughnut) >> 1