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);

  • 相关阅读:
    第十章 Ingress
    第九章 Service
    第八章 资源控制器
    第一章 Xshell5评估期已过问题
    第七章 yaml格式
    第六章 资源清单
    第五章 配置私有仓库Harbor
    第四章 K8s部署安装
    36 SpringBoot 在系统配置文件中动态加载配置
    Java 集合、数组 任意个数数字相加等于一个指定的数
  • 原文地址:https://www.cnblogs.com/fellow1988/p/10004541.html
Copyright © 2011-2022 走看看