Module: Convolver

Defined in:
lib/convolver.rb,
lib/convolver/version.rb,
lib/convolver/operation_plan.rb,
lib/convolver/operation_shapes.rb,
lib/convolver/signal_extension.rb,
lib/convolver/operation_options.rb,
ext/convolver/convolver.c

Overview

Internal planning and extension support for Convolver's public operations.

Constant Summary collapse

MAX_RANK =

Maximum number of dimensions supported by the implementations.

16
VERSION =

Current gem version.

'2.0.0'

Class Method Summary collapse

Class Method Details

.convolve(signal, kernel, mode: :valid, boundary: :constant, fill_value: UNSPECIFIED_FILL, origin: 0) ⇒ Numo::SFloat

Chooses the likely fastest cross-correlation implementation.

Parameters:

  • signal (Numo::NArray) —

    input values

  • kernel (Numo::NArray) —

    correlation kernel

  • mode (:valid, :same, :full) (defaults to: :valid) —

    returned output extent

  • boundary (:constant, :nearest, :reflect, :mirror, :wrap) (defaults to: :constant) —

    signal extension

  • fill_value (Numeric) (defaults to: UNSPECIFIED_FILL) —

    constant extension value

  • origin (Integer, Array<Integer>) (defaults to: 0) —

    kernel origin shift

Returns:

  • (Numo::SFloat) —

    cross-correlation result

Raises:

  • (ArgumentError) —

    if inputs or options are incompatible



26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
 
# File 'lib/convolver.rb', line 26

def convolve(signal, kernel, mode: :valid, boundary: :constant,
             fill_value: UNSPECIFIED_FILL, origin: 0)
  plan = operation_plan(signal, kernel, mode:, boundary:, fill_value:, origin:)
  options = operation_options(mode, boundary, fill_value, origin)

  return convolve_basic(signal, kernel, **options) if plan.extended_size < 1000 || kernel.size < 100

  basic_time_predicted = predict_convolve_basic_time(signal, kernel, **options)
  return convolve_basic(signal, kernel, **options) if basic_time_predicted < 0.1

  fft_time_predicted = predict_convolve_fft_time(signal, kernel, **options)
  return convolve_fft(signal, kernel, **options) if fft_time_predicted < 2 * basic_time_predicted

  convolve_basic(signal, kernel, **options)
end

.convolve_basic(signal, kernel, mode: :valid, boundary: :constant, fill_value: UNSPECIFIED_FILL, origin: 0) ⇒ Numo::SFloat

Uses the direct native valid primitive after applying the requested signal extension in Ruby.

Returns:

  • (Numo::SFloat) —

    cross-correlation result



46
47
48
49
50
 
# File 'lib/convolver.rb', line 46

def convolve_basic(signal, kernel, mode: :valid, boundary: :constant,
                   fill_value: UNSPECIFIED_FILL, origin: 0)
  plan = operation_plan(signal, kernel, mode:, boundary:, fill_value:, origin:)
  convolve_basic_valid(plan.extend_signal(signal), kernel)
end

.convolve_basic_valid(signal, kernel) ⇒ Numo::SFloat

Calculates a valid cross-correlation using the direct native implementation.

Returns valid cross-correlation result.

Parameters:

  • signal (Numo::NArray) —

    input values

  • kernel (Numo::NArray) —

    correlation kernel

Returns:

  • (Numo::SFloat) —

    valid cross-correlation result



25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
 
# File 'ext/convolver/convolver.c', line 25

static VALUE convolver_convolve_basic_valid(VALUE self, VALUE signal, VALUE kernel) {
  volatile VALUE signal_value;
  volatile VALUE kernel_value;
  volatile VALUE result_value;
  narray_t *signal_narray;
  narray_t *kernel_narray;
  int rank;
  int i;
  size_t signal_shape[LARGEST_RANK];
  size_t kernel_shape[LARGEST_RANK];
  size_t result_shape[LARGEST_RANK];
  size_t numo_result_shape[LARGEST_RANK];

  (void)self;

  if (!rb_obj_is_kind_of(signal, numo_cNArray) || !rb_obj_is_kind_of(kernel, numo_cNArray)) {
    rb_raise(rb_eArgError, "signal and kernel must be Numo::NArray values");
  }

  signal_value = rb_funcall(numo_cSFloat, rb_intern("cast"), 1, signal);
  kernel_value = rb_funcall(numo_cSFloat, rb_intern("cast"), 1, kernel);
  if (!RTEST(na_check_contiguous(signal_value))) {
    signal_value = na_copy(signal_value);
  }
  if (!RTEST(na_check_contiguous(kernel_value))) {
    kernel_value = na_copy(kernel_value);
  }
  GetNArray(signal_value, signal_narray);
  GetNArray(kernel_value, kernel_narray);

  if (signal_narray->size == 0 || kernel_narray->size == 0) {
    rb_raise(rb_eArgError, "signal and kernel must not be empty");
  }

  if (signal_narray->ndim != kernel_narray->ndim) {
    rb_raise(rb_eArgError, "signal and kernel must have equal rank");
  }
  if (signal_narray->ndim > LARGEST_RANK) {
    rb_raise(rb_eArgError, "maximum supported rank is %d", LARGEST_RANK);
  }

  rank = signal_narray->ndim;
  copy_shape(rank, signal_narray->shape, signal_shape);
  copy_shape(rank, kernel_narray->shape, kernel_shape);

  for (i = 0; i < rank; i++) {
    if (signal_shape[i] < kernel_shape[i]) {
      rb_raise(rb_eArgError, "kernel must not be larger than signal in any dimension");
    }
    result_shape[i] = signal_shape[i] - kernel_shape[i] + 1;
    numo_result_shape[rank - i - 1] = (size_t)result_shape[i];
  }

  result_value = nary_new(numo_cSFloat, rank, numo_result_shape);
  convolve_raw(
    rank, signal_shape, (float *)na_get_pointer_for_read(signal_value),
    rank, kernel_shape, (float *)na_get_pointer_for_read(kernel_value),
    rank, result_shape, (float *)na_get_pointer_for_write(result_value)
  );

  return result_value;
}

.convolve_fft(signal, kernel, mode: :valid, boundary: :constant, fill_value: UNSPECIFIED_FILL, origin: 0) ⇒ Numo::SFloat

Uses PocketFFT to calculate the requested cross-correlation.

Periodic same-sized results use a circular transform. Other combinations use the shared extension plan and PocketFFT's linear convolution.

Returns:

  • (Numo::SFloat) —

    cross-correlation result



58
59
60
61
62
63
64
 
# File 'lib/convolver.rb', line 58

def convolve_fft(signal, kernel, mode: :valid, boundary: :constant,
                 fill_value: UNSPECIFIED_FILL, origin: 0)
  plan = operation_plan(signal, kernel, mode:, boundary:, fill_value:, origin:)
  return convolve_fft_wrap(signal, kernel, plan) if plan.wrap?

  convolve_fft_valid(plan.extend_signal(signal), kernel)
end

.convolve_fftw3(signal, kernel, mode: :valid, boundary: :constant, fill_value: UNSPECIFIED_FILL, origin: 0) ⇒ Numo::SFloat

Deprecated.

Use convolve_fft; Convolver no longer uses FFTW3.

Compatibility alias for the former FFTW3-backed implementation.

Returns:

  • (Numo::SFloat) —

    cross-correlation result



70
71
72
73
74
75
 
# File 'lib/convolver.rb', line 70

def convolve_fftw3(signal, kernel, mode: :valid, boundary: :constant,
                   fill_value: UNSPECIFIED_FILL, origin: 0)
  warn 'Convolver.convolve_fftw3 is deprecated; use .convolve_fft instead', uplevel: 1
  options = operation_options(mode, boundary, fill_value, origin)
  convolve_fft(signal, kernel, **options)
end

.predict_convolve_basic_time(signal, kernel, mode: :valid, boundary: :constant, fill_value: UNSPECIFIED_FILL, origin: 0) ⇒ Float

Estimates the relative cost of convolve_basic for the requested options.

Returns:

  • (Float) —

    machine-specific relative cost estimate



91
92
93
94
95
96
97
 
# File 'lib/convolver.rb', line 91

def predict_convolve_basic_time(signal, kernel, mode: :valid, boundary: :constant,
                                fill_value: UNSPECIFIED_FILL, origin: 0)
  plan = operation_plan(signal, kernel, mode:, boundary:, fill_value:, origin:)
  operations = plan.result_size * plan.extended_size * kernel.size
  operations += plan.extended_size unless plan.valid?
  4.54e-12 * operations
end

.predict_convolve_fft_time(signal, kernel, mode: :valid, boundary: :constant, fill_value: UNSPECIFIED_FILL, origin: 0) ⇒ Float

Estimates the relative cost of convolve_fft for the requested options.

Returns:

  • (Float) —

    machine-specific relative cost estimate



80
81
82
83
84
85
86
 
# File 'lib/convolver.rb', line 80

def predict_convolve_fft_time(signal, kernel, mode: :valid, boundary: :constant,
                              fill_value: UNSPECIFIED_FILL, origin: 0)
  plan = operation_plan(signal, kernel, mode:, boundary:, fill_value:, origin:)
  transform_size = plan.wrap? ? signal.size : plan.linear_fft_size(kernel.shape)
  transform_cost = 16 * 4.55e-08 * transform_size * Math.log(transform_size)
  transform_cost + (4.55e-08 * fft_preparation_size(plan, signal, kernel))
end