zoukankan      html  css  js  c++  java
  • soundtouch change rate matlab implementation

    soundtouch implement of changing rate in a way same with resample(SRC).

     When rate < 1, it need interpolate sample. and delete samples when rate > 1.

    After  interpolation, there may be introduce high frequency. To avoid aliase, we usally apply a low pass filter on interpolated signal.

    For the case of deleting samples, some signal may contains high frequecy, the difference between sample may be very sharp. Therefore, we use low pass filter before deleting samples.

    When design fir low pass filter, we usaully use "w" (rad/s) as the main parameter instead freqence "f".

    The relationship of "w" and "f" as following:

    w = 2*pi * f/ fs;

    For cutoff frequency, wc = 2*pi * fc / fs = 2 * pi * fc/(2*fn) =  pi * fc / fn, where fc is cutoff frequency, fn is Nywuist frequency.

    For fc = = fn, wc = pi;

    We use cubic method to interplate sample, the principle of cubic interpolation refer to https://www.paulinternet.nl/?page=bicubic

    %calc low pass filter coefficient. The low pass filter based on sinc function with hamming window.

    function coeff = calCoeffs(cutoffFreq, len)

    coeff = zeros(len ,1);

    wc = 2 * pi * cutoffFreq;

    tempCoeff = 2 * pi / len;

    sum = 0;

    for i = 0 : 1 : len -1

      cntTemp = (i - len/2);

      temp = cntTemp  * wc;

      % sinc function

      if temp ~=0

        h = sin(temp) / temp;

      else

        h = 1;

      end

      %hamming window

      w = 0.54 + 0.46 * cos(tempCoeff * cntTemp);

      coeff(i+1) = w * h;

      sum = sum + coeff(i+1);

    end

    coeff = coeff / sum;

    end

    function output = firfilter(input, coeff)

    inputLen = length(input(:, 1));

    filterLen = length(coeff(:, 1));

    output = zeros(inputLen ,1 );

    outputLen = inputLen - filterLen;

    for i = 1: 1: outputLen

      inpos = i;

      sum = 0;

      for j = 1:1:filterLen

        sum = sum + input(inpos ,1) * coeff(j, 1);

        inpos = inpos + 1;

      end

      output(i, 1) = sum;

    end

    end

    function output = cubicInterpolation(input, rate)

    inputLen = length(input(:,1));

    outputLen = floor(inputLen / rate);

    output = zeros(outputLen ,1);

    inputIdx = 1;

    fract = 0;

    outputIdx = 1;

    while inputIdx < inputLen - 4

      x1 = fract;

      x2 = x1 * x1;

      x3 = x1 * x2;

      p0 = input(inputIdx , 1);

      p1 = input(inputIdx + 1 , 1);

         p2 = input(inputIdx + 2, 1);

      p3 = input(inputIdx + 3, 1);

      output(outputIdx ,1) = (-0.5*p0 + 1.5*p1 -1.5 *p2 + 0.5*p3) * x3 +(p0 - 2.5*p1 + 2*p2 -0.5*p3) *x2 + (-0.5*p0 + 0.5*p2) * x1 + p1;

      outputIdx = outputIdx + 1;

      fract = fract + rate;

      whole = floor(fract);

      fract = fract - whole;

      inputIdx = inputIdx + whole;

     end

    end

    function output = linearInterpolation(input, rate)

    inputLen = length(input(:,1));

    outputLen = floor(inputLen / rate);

    output = zeros(outputLen ,1);

    inputIdx = 1;

    fract = 0;

    outputIdx = 1;

    while inputIdx < inputLen - 4

      p0 = input(inputIdx , 1);

      p1 = input(inputIdx + 1 , 1);

      output(outputIdx ,1) = (1-fract) * po + fract * p1;

      outputIdx = outputIdx + 1;

      fract = fract + rate;

      whole = floor(fract);

      fract = fract - whole;

      inputIdx = inputIdx + whole;

     end

    end

    function output = changeRate(input, rate, interpMethod)

    inputLen = length(input(:, 1));

    outputLen = floor(inputLen / rate);

    output = zeros(outputLen, 1);

    if rate > 1

      cutoffFreq = 0.5 / rate;

    else

      cutoffFreq = 0.5 * rate;

    end

    filterLen = 64;

    coeff = calCoeffs(cutoffFreq, filterLen);

    if rate < 1

      %slow down, need interpolation first;

      if strcmp(interMethod, 'cubic')

        output = cubicInterpolation(input, rate);

      else

        output = linearInterpolation(input, rate);

       end

      output = firfilter(output, coeff);

    else

      %fast, need filter out the high freqency, then delete samples

      output = firfilter(input, coeff);

      if strcmp(interMethod, 'cubic')

        output = cubicInterpolation(output, rate);

      else

        output = linearInterpolation(output, rate);

       end

    end

    end

    main.m:

    clc;

    clear all;

    [input fs] = wavread('input.wav');

    %if do SRC, rate = inputfs / outputfs;

    rate = 0.5;

    output = changeRate(input, rate, 'cubic');

    wavwrite(output, fs, 'output.wav);

  • 相关阅读:
    数学--数论---欧拉筛 模板
    数学--数论--HDU 5019 revenge of GCD
    数学--数论--HDU 2674 沙雕题
    数学--数论--HDU 2104 丢手绢(离散数学 mod N+ 剩余类 生成元)+(最大公约数)
    数学--数论-- HDU -- 2854 Central Meridian Number (暴力打表)
    数学--数论--HDU 2197 本原串 (推规律)
    数学--数论--HDU 2582 F(N) 暴力打表找规律
    C# 命名规范
    文本分类学习 (十)构造机器学习Libsvm 的C# wrapper(调用c/c++动态链接库)
    文本分类学习 (九)SVM入门之拉格朗日和KKT条件
  • 原文地址:https://www.cnblogs.com/fellow1988/p/10004541.html
Copyright © 2011-2022 走看看