purify::kernels Namespace Reference

Enumerations
enum class	kernel {   kb , gauss , box , pswf ,   kbmin , gauss_alt , kb_presample }

Functions
t_real	kaiser_bessel (const t_real x, const t_real J)
	Kaiser-Bessel kernel. More...
t_real	kaiser_bessel_general (const t_real x, const t_real J, const t_real alpha)
	More general Kaiser-Bessel kernel. More...
t_real	ft_kaiser_bessel_general (const t_real x, const t_real J, const t_real alpha)
	Fourier transform of more general Kaiser-Bessel kernel. More...
t_real	ft_kaiser_bessel (const t_real x, const t_real J)
	Fourier transform of kaiser bessel kernel. More...
t_real	gaussian (const t_real x, const t_real J)
	Gaussian kernel. More...
t_real	ft_gaussian (const t_real x, const t_real J)
	Fourier transform of Gaussian kernel. More...
t_real	calc_for_pswf (const t_real x, const t_real J, const t_real alpha)
	Calculates Horner's Rule the standard PSWF for radio astronomy, with a support of J = 6 and alpha = 1. More...
t_real	pswf (const t_real x, const t_real J)
	PSWF kernel. More...
t_real	ft_pswf (const t_real x, const t_real J)
	Fourier transform of PSWF kernel. More...
std::vector< t_real >	kernel_samples (const t_int total_samples, const std::function< t_real(t_real)> kernelu)
	Calculates samples of a kernel. More...
t_real	kernel_zero_interp (const std::vector< t_real > &samples, const t_real x, const t_real J)
	zeroth order interpolates from samples of kernel More...
t_real	kernel_linear_interp (const Vector< t_real > &samples, const t_real x, const t_real J)
	linearly interpolates from samples of kernel More...
t_real	pill_box (const t_real x, const t_real J)
	Box car function for kernel. More...
t_real	ft_pill_box (const t_real x, const t_real J)
	Fourier transform of box car function, a Sinc function. More...
t_real	gaussian_general (const t_real x, const t_real J, const t_real sigma)
	Fourier transform of general Gaussian kernel. More...
t_real	ft_gaussian_general (const t_real x, const t_real J, const t_real sigma)
	Fourier transform of general Gaussian kernel. More...

Variables
const std::map< std::string, kernel >	kernel_from_string

Enumeration Type Documentation

kernel

enum purify::kernels::kernel

strong

Enumerator
kb
gauss
box
pswf
kbmin
gauss_alt
kb_presample

Definition at line 15 of file kernels.h.

{ kb, gauss, box, pswf, kbmin, gauss_alt, kb_presample };

Function Documentation

calc_for_pswf()

t_real purify::kernels::calc_for_pswf	(	const t_real	x,
		const t_real	J,
		const t_real	alpha
	)

Calculates Horner's Rule the standard PSWF for radio astronomy, with a support of J = 6 and alpha = 1.

Parameters

[in]	eta0	value to evaluate
[in]	J	support size of gridding kernel
[in]	alpha	type of special PSWF to calculate

The tailored prolate spheroidal wave functions for gridding radio astronomy. Details are explained in Optimal Gridding of Visibility Data in Radio Astronomy, F. R. Schwab 1983.

Definition at line 71 of file kernels.cc.

                                                                            {
  // polynomial coefficients for prolate spheriodal wave function rational approximation
  const std::array<t_real, 6> p1 = {
      {8.203343e-2, -3.644705e-1, 6.278660e-1, -5.335581e-1, 2.312756e-1, 2 * 0.0}};
  const std::array<t_real, 6> p2 = {
      {4.028559e-3, -3.697768e-2, 1.021332e-1, -1.201436e-1, 6.412774e-2, 2 * 0.0}};
  const std::array<t_real, 3> q1 = {{1., 8.212018e-1, 2.078043e-1}};
  const std::array<t_real, 3> q2 = {{1., 9.599102e-1, 2.918724e-1}};
 
  if (J != 6 or alpha != 1) {
    return 0;
  }
  // Calculating numerator and denominator using Horner's rule.
  // PSWF = numerator / denominator
 
  auto fraction = [](t_real eta, std::array<t_real, 6> const &p, std::array<t_real, 3> const &q) {
    auto const p_size = sizeof(p) / sizeof(p[0]) - 1;
    auto const q_size = sizeof(q) / sizeof(q[0]) - 1;
 
    auto numerator = p[p_size];
    for (auto i = decltype(p_size){1}; i <= p_size; ++i)
      numerator = eta * numerator + p[p_size - i];
 
    auto denominator = q[q_size];
    for (auto i = decltype(q_size){1}; i <= q_size; ++i)
      denominator = eta * denominator + q[q_size - i];
 
    return numerator / denominator;
  };
  if (0 <= std::abs(eta0) and std::abs(eta0) <= 0.75)
    return fraction(eta0 * eta0 - 0.75 * 0.75, p1, q1);
  if (0.75 < std::abs(eta0) and std::abs(eta0) <= 1) return fraction(eta0 * eta0 - 1 * 1, p2, q2);
 
  return 0;
}

Referenced by ft_pswf(), and pswf().

ft_gaussian()

t_real purify::kernels::ft_gaussian	(	const t_real	x,
		const t_real	J
	)

Fourier transform of Gaussian kernel.

Definition at line 60 of file kernels.cc.

                                                   {
  /*
     Fourier transform of guassian gridding kernel
 
     x:: value to evaluate
     J:: support size of gridding kernel
     */
  t_real sigma = 0.31 * std::pow(J, 0.52);
  return ft_gaussian_general(x, J, sigma);
}

References ft_gaussian_general().

Referenced by purify::create_kernels(), and purify::create_radial_ftkernel().

ft_gaussian_general()

t_real purify::kernels::ft_gaussian_general	(	const t_real	x,
		const t_real	J,
		const t_real	sigma
	)

Fourier transform of general Gaussian kernel.

Definition at line 233 of file kernels.cc.

                                                                               {
  /*
     Fourier transform of Gaussian gridding kernel
 
     x:: value to evaluate
     J:: support size of gridding kernel
     sigma:: standard deviation of Gaussian kernel (in pixels)
     */
 
  t_real a = x * sigma;
  return std::exp(a * a * 2);
}

Referenced by purify::create_kernels(), and ft_gaussian().

ft_kaiser_bessel()

t_real purify::kernels::ft_kaiser_bessel	(	const t_real	x,
		const t_real	J
	)

Fourier transform of kaiser bessel kernel.

Definition at line 40 of file kernels.cc.

                                                        {
  /*
     Fourier transform of kaiser bessel gridding kernel
     */
 
  t_real alpha = 2.34 * J;  // value said to be optimal in Fessler et. al. 2003
  return ft_kaiser_bessel_general(x, J, alpha);
}

References ft_kaiser_bessel_general().

Referenced by purify::create_kernels(), purify::create_radial_ftkernel(), main(), and TEST_CASE().

ft_kaiser_bessel_general()

t_real purify::kernels::ft_kaiser_bessel_general	(	const t_real	x,
		const t_real	J,
		const t_real	alpha
	)

Fourier transform of more general Kaiser-Bessel kernel.

Definition at line 25 of file kernels.cc.

                                                                                    {
  /*
     Fourier transform of kaiser bessel gridding kernel
     */
 
  t_complex eta = std::sqrt(
      static_cast<t_complex>((constant::pi * x * J) * (constant::pi * x * J) - alpha * alpha));
  const t_real normalisation = J / boost::math::cyl_bessel_i(0, alpha);
 
  return std::real(std::sin(eta) / eta) *
         normalisation;  // simple way of doing the calculation, the
  // boost bessel funtions do not support
  // complex valued arguments
}

References purify::constant::pi.

Referenced by purify::create_kernels(), purify::create_radial_ftkernel(), and ft_kaiser_bessel().

ft_pill_box()

t_real purify::kernels::ft_pill_box	(	const t_real	x,
		const t_real	J
	)

Fourier transform of box car function, a Sinc function.

Definition at line 211 of file kernels.cc.

                                                   {
  /*
     Gaussian gridding kernel
 
     x:: value to evaluate
     J:: support size
     */
  return boost::math::sinc_pi(J * x * constant::pi);
}

References purify::constant::pi.

Referenced by purify::create_kernels(), and purify::create_radial_ftkernel().

ft_pswf()

t_real purify::kernels::ft_pswf	(	const t_real	x,
		const t_real	J
	)

Fourier transform of PSWF kernel.

Definition at line 124 of file kernels.cc.

                                               {
  /*
     Calculates the fourier transform of the standard PSWF for radio astronomy, with a support of J
     = 6 and alpha = 1.
 
     x:: value to evaluate
     J:: support size of gridding kernel
 
     The tailored prolate spheroidal wave functions for gridding radio astronomy.
     Details are explained in Optimal Gridding of Visibility Data in Radio
     Astronomy, F. R. Schwab 1983.
 
*/
 
  const t_real alpha = 1;
  const t_real eta0 = 2 * x;
 
  return calc_for_pswf(eta0, J, alpha) * std::sqrt(5.343);
}

References calc_for_pswf().

Referenced by purify::create_kernels(), and purify::create_radial_ftkernel().

gaussian()

t_real purify::kernels::gaussian	(	const t_real	x,
		const t_real	J
	)

Gaussian kernel.

Definition at line 49 of file kernels.cc.

                                                {
  /*
     Guassian gridding kernel
 
     x:: value to evaluate
     J:: support size
     */
  t_real sigma = 0.31 * std::pow(J, 0.52);  // Optimal sigma according to fessler et al.
  return gaussian_general(x, J, sigma);
}

References gaussian_general().

Referenced by purify::create_kernels(), and purify::create_radial_ftkernel().

gaussian_general()

t_real purify::kernels::gaussian_general	(	const t_real	x,
		const t_real	J,
		const t_real	sigma
	)

Fourier transform of general Gaussian kernel.

Definition at line 220 of file kernels.cc.

                                                                            {
  /*
     Gaussian gridding kernel
 
     x:: value to evaluate
     J:: support size
     sigma:: standard deviation of Gaussian kernel (in pixels)
     */
 
  t_real a = x / sigma;
  return std::exp(-a * a * 0.5) / (sigma * std::sqrt(2 * constant::pi));
}

References purify::constant::pi.

Referenced by purify::create_kernels(), and gaussian().

kaiser_bessel()

t_real purify::kernels::kaiser_bessel	(	const t_real	x,
		const t_real	J
	)

Kaiser-Bessel kernel.

Definition at line 7 of file kernels.cc.

                                                     {
  /*
     kaiser bessel gridding kernel
     */
  t_real alpha = 2.34 * J;  // value said to be optimal in Fessler et. al. 2003
  return kaiser_bessel_general(x, J, alpha);
}

References kaiser_bessel_general().

Referenced by purify::create_kernels(), purify::create_radial_ftkernel(), and TEST_CASE().

kaiser_bessel_general()

t_real purify::kernels::kaiser_bessel_general	(	const t_real	x,
		const t_real	J,
		const t_real	alpha
	)

More general Kaiser-Bessel kernel.

Definition at line 15 of file kernels.cc.

                                                                                 {
  /*
     kaiser bessel gridding kernel
     */
  if (2 * x > J) return 0;
  t_real a = 2 * x / J;
  return boost::math::cyl_bessel_i(0, std::real(alpha * std::sqrt(1 - a * a))) /
         boost::math::cyl_bessel_i(0, alpha);
}

Referenced by purify::create_kernels(), purify::create_radial_ftkernel(), and kaiser_bessel().

kernel_linear_interp()

t_real purify::kernels::kernel_linear_interp	(	const Vector< t_real > &	samples,
		const t_real	x,
		const t_real	J
	)

linearly interpolates from samples of kernel

Definition at line 171 of file kernels.cc.

                                                                                           {
  /*
     Calculates kernel using linear interpolation between pre-calculated samples. (see Rapid
     gridding reconstruction with a minimal oversampling ratio, Beatty et. al. 2005)
     */
  t_int total_samples = samples.size();
 
  t_real i_effective = (x + J / 2) * total_samples / J;
 
  t_int i_0 = floor(i_effective);
  t_int i_1 = ceil(i_effective);
  // case where i_effective is a sample point
  if (std::abs(i_0 - i_1) == 0) {
    return samples(i_0);
  }
  // linearly interpolate from nearest neighbour
  t_real y_0;
  t_real y_1;
  if (i_0 < 0 or i_0 >= total_samples) {
    y_0 = 0;
  } else {
    y_0 = samples(i_0);
  }
  if (i_1 < 0 or i_1 >= total_samples) {
    y_1 = 0;
  } else {
    y_1 = samples(i_1);
  }
  t_real output = y_0 + (y_1 - y_0) / (i_1 - i_0) * (i_effective - i_0);
  return output;
}

kernel_samples()

std::vector< t_real > purify::kernels::kernel_samples	(	const t_int	total_samples,
		const std::function< t_real(t_real)>	kernelu
	)

Calculates samples of a kernel.

Definition at line 144 of file kernels.cc.

                                                                              {
  /*
     Pre-calculates samples of a kernel, that can be used with linear interpolation (see Rapid
     gridding reconstruction with a minimal oversampling ratio, Beatty et. al. 2005)
     */
  std::vector<t_real> samples(total_samples);
  for (Vector<t_real>::Index i = 0; i < static_cast<t_int>(total_samples); ++i) {
    samples[i] = kernelu(static_cast<t_real>(i) / total_samples);
  }
  return samples;
}

Referenced by purify::create_kernels(), and purify::operators::init_on_the_fly_gridding_matrix_2d().

kernel_zero_interp()

t_real purify::kernels::kernel_zero_interp	(	const std::vector< t_real > &	samples,
		const t_real	x,
		const t_real	J
	)

zeroth order interpolates from samples of kernel

Definition at line 156 of file kernels.cc.

                                                                                            {
  /*
     Calculates kernel using linear interpolation between pre-calculated samples. (see Rapid
     gridding reconstruction with a minimal oversampling ratio, Beatty et. al. 2005)
     */
  const t_int total_samples = samples.size();
 
  const t_real i_effective = 2 * std::abs(x) * total_samples / J;
 
  const t_real i_0 = std::floor(i_effective);
  if (i_0 < 0 or i_0 >= total_samples) return 0.;
 
  return samples[static_cast<t_int>(i_0)];
}

Referenced by purify::create_kernels().

pill_box()

t_real purify::kernels::pill_box	(	const t_real	x,
		const t_real	J
	)

Box car function for kernel.

Definition at line 202 of file kernels.cc.

                                                {
  /*
     Gaussian gridding kernel
 
     x:: value to evaluate
     J:: support size
     */
  return 1. / static_cast<t_real>(J);
}

Referenced by purify::create_kernels(), and purify::create_radial_ftkernel().

pswf()

t_real purify::kernels::pswf	(	const t_real	x,
		const t_real	J
	)

PSWF kernel.

Definition at line 107 of file kernels.cc.

                                            {
  /*
     Calculates the standard PSWF for radio astronomy, with a support of J = 6 and alpha = 1.
 
     x:: value to evaluate
     J:: support size of gridding kernel
 
     The tailored prolate spheroidal wave functions for gridding radio astronomy.
     Details are explained in Optimal Gridding of Visibility Data in Radio
     Astronomy, F. R. Schwab 1983.
 
*/
  const t_real eta0 = 2 * x / J;
  const t_real alpha = 1;
  return calc_for_pswf(eta0, J, alpha) * std::pow(1 - eta0 * eta0, alpha);
}

References calc_for_pswf().

Referenced by purify::create_kernels(), and purify::create_radial_ftkernel().

Variable Documentation

kernel_from_string

const std::map<std::string, kernel> purify::kernels::kernel_from_string

Initial value:

= {{"kb", kernel::kb},
                                                          {"gauss", kernel::gauss},
                                                          {"box", kernel::box},
                                                          {"pswf", kernel::pswf},
                                                          {"kbmin", kernel::kbmin},
                                                          {"kb_presample", kernel::kb_presample},
                                                          {"gauss_alt", kernel::gauss_alt}}

Definition at line 16 of file kernels.h.

Referenced by BENCHMARK_DEFINE_F(), createMeasurementOperator(), getInputData(), main(), padmm(), and TEST_CASE().