Module: Savgol

Included in:
Array
Defined in:
lib/savgol.rb,
lib/savgol/array.rb,
lib/savgol/version.rb

Constant Summary collapse

VERSION = 
"0.2.1"

Instance Method Summary collapse

Instance Method Details

#savgol(window_size, order, deriv = 0, check_args = false) ⇒ Object 



5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
 
# File 'lib/savgol/array.rb', line 5

def savgol(window_size, order, deriv=0, check_args=false)
    
  # check arguments
  if !window_size.is_a?(Integer) || window_size.abs != window_size || window_size % 2 != 1 || window_size < 1
    raise ArgumentError, "window_size size must be a positive odd integer" 
  end
  if !order.is_a?(Integer) || order < 0
    raise ArgumentError, "order must be an integer >= 0"
  end
  if window_size < order + 2
    raise ArgumentError, "window_size is too small for the polynomials order"
  end
  
  half_window = (window_size -1) / 2
  weights = sg_weights(half_window, order, deriv)
  ar = sg_pad_ends(half_window)
  sg_convolve(ar, weights)
end

#sg_convolve(data, weights, mode = :valid) ⇒ Object 



24
25
26
27
28
 
# File 'lib/savgol/array.rb', line 24

def sg_convolve(data, weights, mode=:valid)
  data.each_cons(weights.size).map do |ar|
    ar.zip(weights).map {|pair| pair[0] * pair[1] }.reduce(:+)
  end
end

#sg_pad_ends(half_window) ⇒ Object

pads the ends with the reverse, geometric inverse sequence



31
32
33
34
35
36
37
38
39
40
 
# File 'lib/savgol/array.rb', line 31

def sg_pad_ends(half_window)
  start = self[1..half_window]
  start.reverse!
  start.map! {|v| self[0] - (v - self[0]).abs }

  fin = self[(-half_window-1)...-1]
  fin.reverse!
  fin.map! {|v| self[-1] + (v - self[-1]).abs }
  start.push(*self, *fin)
end

#sg_weights(half_window, order, deriv = 0) ⇒ Object

returns an object that will convolve with the padded array



43
44
45
46
47
48
49
50
51
 
# File 'lib/savgol/array.rb', line 43

def sg_weights(half_window, order, deriv=0)
  # byebug
  mat = Matrix[ *(-half_window..half_window).map {|k| (0..order).map {|i| k**i }} ]
  # Moore-Penrose psuedo-inverse without SVD (not so precize)
  # A' = (A.t * A)^-1 * A.t
  pinv_matrix = Matrix[*(mat.transpose*mat).to_a].inverse * Matrix[*mat.to_a].transpose
  pinv = Matrix[*pinv_matrix.to_a]
  pinv.row(deriv).to_a
end