ディジタル信号処理による信号発生とエイリアス(11)

前回の「gen_saw2.sce」プログラムに、オーバーサンプリング/デシメーションによる波形である「saw4s.wav」の生成を追加しました。

// gen_saw2.sce - Saw wave gen by Scilab
// 
// globals
//
global nperi harm maxlevel n_ovs;
pi = 4.0 * atan(1.0); // constant pi
nharm = 5;   // number of harmonics
ncyc = 1024; // cycles per block
nrep = 48;   // block repeat count
fs = 48000;  // sampling frequency
fre = 1023.0 / ncyc; // normalized freq
nperi = 11;  // waveform period 
n_ovs = 2;   // oversampling factor
harm = (1 : nharm); // harmonics vector
maxlevel = (6.0 / 8.0); // = 0x6000/0x8000
rep = ones(1, nrep); // repeat vector
x = (0 : ncyc-1); // sampling point vector
x = fre * x; // multiply by nomalized freq
// sampled ideal saw wave
function r=saw2(x)
  global nperi maxlevel;
  y = ((x) ./ nperi);
  y = y - floor(y);
  r = 2.0 * maxlevel * (y - 0.5);
endfunction
//
// Band limited saw by Fourier series
function r=saw3(x)
  global nperi harm maxlevel;
  y = ((x) ./ nperi);
  y = y - floor(y);
  y = ((harm.') * y);
  s = ((1.0) ./ harm) * sin(2.0 * pi * y);
  r = (-2.0 * maxlevel / pi) * s;
endfunction
//
// N times oversampled/decimated ideal saw wave
function r=saw4(x)
  global nperi maxlevel n_ovs;
  y = x; // original sampling points
  a = (x(2) - x(1)); // sampling interval
  for k = (1:n_ovs-1); // add ovs offset
    y = [y; (k * a / n_ovs) + x];
  end;
  y = ((y) ./ nperi);  
  y = 2.0 * (y - floor(y) - 0.5);
  y = sum(y, 1) / n_ovs; // comb filter
  r = maxlevel * y;
endfunction
//
// main
//
// Lch = sampled ideal sawtooth wave
// Rch = band limited saw wave by Fourier series
lr = [saw2(x); saw3(x)];
wavwrite(kron(rep, lr), fs, "saw2s.wav");
//
// Lch = N times oversampled/decimated ideal saw
// Rch = band limited saw wave by Fourier series
lr = [saw4(x); saw3(x)];
wavwrite(kron(rep, lr), fs, "saw4s.wav");

コメントの「//」 (スラッシュ、スラッシュ) を「#」 (ポンド記号) に置換し、最後の行の wavwrite() 関数の第一引数を転置すれば、GNU Octave 3.0.0 でも実行できます。