ProFit (version 1.3.3)

profitBenchmark: Benchmark profile integration and image convolution using libprofit.

Description

This function will benchmark integration of surface brightness profiles and/or convolution of an image with a kernel (usually a point spread function) using libprofit. It returns a data frame with available integration/convolution methods, the most efficient method(s), and (optionally) more detailed results including convolution accuracy. It is called by profitSetupData by default and the convolver is used in profitMakeModel.

Usage

profitBenchmark(image, methods=NULL, psf=NULL,
  modellist=NULL, finesample=1L, calcregion=NULL, nbench=1,
  benchconvolution=is.matrix(psf),
  precisions=c("double"), omp_threads=1,
  openclenvs = profitGetOpenCLEnvs(make.envs = TRUE),
  reference = "brute", reusepsffft = TRUE, fft_effort=0,
  returnimages = FALSE)

Arguments

image

A matrix containing the image to benchmark convolution for. It should already be padded by half of the PSF width on either side to ensure that the convolved model can be cropped to the same size as the data. If no image is supplied, the user must supply imagedim.

methods

List of strings specifying which methods to test. Methods must be amongst those returned by profitAvailableConvolvers.

psf

A matrix containing the PSF image to convolve the model image with. If no PSF is supplied, the user must supply psfdim.

modellist

A valid ProFit modellist for profile integration; see profitMakeModel.

finesample

Integer; the factor to oversample the image by. The default of one does no oversampling.

calcregion

A logical matrix specifying regions of the image to avoid computing convolution for. Can make brute force convolution more efficient if it is not sparse.

nbench

Integer; the number of times to benchmark each method. Repeated convolutions can vary in running time for all kinds of reasons, so nbench = 10 is recommended unless using brute force convolution with very large images and/or kernels.

benchconvolution

Logical; whether to benchmark convolution. Requires a PSF.

precisions

Character; the numerical precision(s) to benchmark. Must be one of "single" or "double".

omp_threads

Integer; the number of OpenMP threads to use for integration/convolution.

openclenvs

A data.frame with information on available OpenCL environments, such as that returned from profitGetOpenCLEnvs.

reference

String; the method to use as the reference result for comparing the accuracy of all other methods. This comparison is not done if reference is not contain in methods.

reusepsffft

Logical specifying whether to re-do the PSF FFT every iteration, which would be necessary if one is fitting the PSF.

fft_effort

The effort level to compute the FFTW plan. FFTW plans can take a very long time to set up, so consider carefully before increasing beyond 0 - particularly if your padded image only has a few large prime factors.

returnimages

Logical; whether to return the convolved image for every method.

Value

List containing:

result

The benchmarking results in a data.frame; see profitGetOpenCLEnvs for more information on the format.

images

List of resulting images for each method.

Notes

TBD.

Details

The function is mainly used to determine the most efficient method for convoling the image with the psf. In situations where the psf has much smaller dimensions than image this will pretty much always be Brute force convolution, but when the psf becomes comparable in size to the image then FFTW is usually faster. In the provided example all three are similar speed. Benchmarks are more difficult to predict when using multiple cores and/or devices.

See Also

profitAvailableConvolvers, profitMakeModel, profitSetupData

Examples

# NOT RUN {
model = list(
	sersic = list(
		xcen   = c(180, 60),
		ycen   = c(90, 10),
		mag = c(15, 13),
		re  = c(14, 5),
		nser  = c(3, 10),
		ang  = c(46, 80),
		axrat  = c(0.4, 0.6),
		box = c(0.5,-0.5)
	)
)

psffwhm=3

# Use OpenCL if available
# Makes a list of available OpenCL environments optionally with double precision if all
# devices support it

openclenvs = profitGetOpenCLEnvs(make.envs=TRUE)
nbench=1L

# Try up to 5L if you're adventurous and don't mind waiting up to a minute for
# single-threaded brute

for(finesample in c(1L:3L))
{
  model.image=profitMakeModel(model=model, dim=rep(200,2), finesample=finesample, returnfine=TRUE)$z
  psf=profitMakeGaussianPSF(fwhm=3*finesample,dim=rep(25*finesample + 1 - mod(finesample,2),2))
  
  # Benchmark model integration:
  bench=profitBenchmark(model.image, modellist=model, nbench=nbench, openclenvs=openclenvs,
    methods=profitAvailableIntegrators())

  #Print relevant results
  print(profitBenchmarkResultStripPointers(bench$result)[
    c("name","env_name","version","dev_name",paste0("tinms.mean_",c("single","double")))])
  
  # Benchmark convolution:
  bench=profitBenchmark(model.image, psf=psf, nbench=nbench, openclenvs=openclenvs,
    methods=profitAvailableConvolvers())
  
  #Print relevant results
  print(profitBenchmarkResultStripPointers(bench$result)[
    c("name","env_name","version","dev_name",paste0("tinms.mean_",c("single","double")))])

  # The old benchmarking method, for reference
  temp=profitBenchmarkConv(model.image, psf = psf, nbench=nbench)
}
# }