win_06.cc

这是一个从音频信号里提取特征参量的程序

// file: $isip/class/algo/Window/win_06.cc
// version: $Id: win_06.cc,v 1.12 2002/05/31 21:57:33 picone Exp $
//

// isip include files
//
#include "Window.h"
#include "Bessel.h"
#include "Chebyshev.h"

// method: setConstants
//
// arguments: 
//  const VectorFloat& values: (input) vector of constants
//
// return: a boolean value indicating status
//
// this method sets the constants of the window
//
boolean Window::setConstants(const VectorFloat& values_a) {
  
  // make sure the correct number of coefficients are specified
  // for the given algorithm
  //
  if (algorithm_d == RECTANGULAR) {
    if (values_a.length() != DEF_RECT_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == BARTLETT) {
    if (values_a.length() != DEF_BART_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == BLACKMAN) {
    if (values_a.length() != DEF_BLCK_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == DOLPH_CHEBYSHEV) {
    if (values_a.length() != DEF_DCHB_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == GAUSSIAN) {
    if (values_a.length() != DEF_GAUS_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == HAMMING) {
    if (values_a.length() != DEF_HAMM_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == HANNING) {
    if (values_a.length() != DEF_HANN_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == KAISER) {
    if (values_a.length() != DEF_KAIS_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == LIFTER) {
    if (values_a.length() != DEF_LIFT_CONSTANTS.length()) {
      return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
    }
  }
  else if (algorithm_d == CUSTOM) {
    return Error::handle(name(), L"setConstants", ERR, __FILE__, __LINE__);
  }
  else {
    return Error::handle(name(), L"setConstants", Error::ENUM,
			 __FILE__, __LINE__);
  }
  
  // set the constants
  //
  constants_d.assign(values_a);
  
  // make the window invalid
  //
  is_valid_d = false;
  
  // exit gracefully
  //
  return true;
}

// method: generateRectangular
//
// arguments: 
//  VectorFloat& output: (output) the window function
//  const VectorFloat& params: (input) window function parameters
//  long num_points: (input) number of points in the window
//
// return: a boolean value indicating status
//
// this method generates a rectangular window function:
//
//  w(n) = { G   n <= |M|
//         { 0   elsewhere
//
// there is one parameter allowed for this window:
//
//  params(0): G, the gain of the window
//
boolean Window::generateRectangular(VectorFloat& output_a,
				    const VectorFloat& params_a,
				    long num_points_a) const {
  
  // for even or odd number of points, do the same thing
  //
  return output_a.assign(num_points_a, params_a(0));
}

// method: generateBartlett
//
// arguments: 
//  VectorFloat& output: (output) the window function
//  const VectorFloat& params: (input) window function parameters
//  long num_points: (input) number of points in the window
//
// return: a boolean value indicating status
//
// this method generates a Bartlett (triangular) window function:
//
//  w(n) = { G (1 - |n| / (M + 1))  n <= |M|
//         { 0                      elsewhere
//
// where M = numn_points_a / 2. a good reference for this window
// can be found at:
//
//  R.A. Roberts and C.T. Mullis, Digital Signal Processing,
//  Addison Wesley Publishing Company, Reading, Massachusetts, USA,
//  1987, pp. 183-188.
//
// there is one parameter allowed for this window:
//
//  params(0): G, the gain of the window
//
// for an odd number of points, the window is defined as above.
// for an even number of points, an interpolation technique is used.
//
boolean Window::generateBartlett(VectorFloat& output_a,
				 const VectorFloat& params_a,
				 long num_points_a) const {
  
  // compute some constants related to the length of the window
  //
  long N = num_points_a;
  long Nm1 = N - 1;
  long Nm2 = Nm1 - 1;
  float scale = (float)Nm2 / (float)Nm1;
  
  // compute the half-window length
  //
  long M = (N - 1) / 2;
  long Mp1 = M + 1;
  
  // create output space
  //
  output_a.setLength(N);
  
  // case 1: an odd number of points
  //
  if ((N % 2) != 0) {
    for (long n = -M; n <= M; n++) {
      output_a(n + M) = (float)params_a(0) * (1 - Integral::abs(n) / (float)(Mp1));
    }
  }
  
  // case 2: an even number of points
  //  resample an (N-1) point window
  //
  else {
    
    // generate a window for (N-1) points
    // sample it at a new time scale 
    //
    for (long n = 0; n < N; n++) {
      float scaled_index = -M + scale * n;
      output_a(n) = (float)params_a(0) *
	(1 - Integral::abs(scaled_index) / (float)(Mp1));
    }
  }    
  
  // exit gracefully
  //
  return true;
}

// method: generateBlackman
//
// arguments: 
//  VectorFloat& output: (output) the window function
//  const VectorFloat& params: (input) window function parameters
//  long num_points: (input) number of points in the window
//
// return: a boolean value indicating status
//
// this method generates a Blackman window function:
//
//  w(n) = { G (a1 + a2 * cos(2*PI*n / (2*M+1))
//         {       + (0.5 - a1) * cos(4*PI*n / (2*M+1))    n <= |M|
//         { 0                                             elsewhere
//
// a good reference for this window can be found at:
//
//  S.K. Mitra, Digital Signal Processing,
//  McGraw-Hill, Boston, Massachusetts, USA, 2001, pp. 452.
//
// there are three parameters allowed for this window:
//
//  params(0): G, the gain of the window
//  params(1): a1, the DC offset
//  params(2): a2, the first harmonic weight
//
// for an odd number of points, the window is defined as above.
// for an even number of points, an interpolation technique is used.
//
boolean Window::generateBlackman(VectorFloat& output_a,
				 const VectorFloat& params_a,
				 long num_points_a) const {
  
  // compute some constants related to the length of the window
  //
  long N = num_points_a;
  long Nm1 = N - 1;
  long Nm2 = Nm1 - 1;
  float scale = (float)Nm2 / (float)Nm1;
  
  // compute the half-window length
  //
  long M = (N - 1) / 2;
  long Mp1 = M + 1;
  
  // create output space
  //
  output_a.setLength(N);
  
  // case 1: an odd number of points
  //
  if ((N % 2) != 0) {
    for (long n = -M; n <= M; n++) {
      float arg1 = 2.0 * Integral::PI * (float)(n) / (float)Mp1;
      float arg2 = 2.0 * arg1;
      output_a(n + M) = (float)params_a(0) *
	((float)params_a(1) + (float)params_a(2) * Integral::cos(arg1) +
	 ((float)params_a(2) - (float)params_a(1)) * Integral::cos(arg2));
    }
  }
  
  // case 2: an even number of points
  //  resample an (N-1) point window
  //
  else {
    
    // generate a window for (N-1) points
    // sample it at a new time scale 
    //
    for (long n = 0; n < N; n++) {
      float scaled_index = -M + scale * n;
      float arg1 = 2.0 * Integral::PI * scaled_index / (float)Mp1;
      float arg2 = 2.0 * arg1;
      output_a(n) = (float)params_a(0) *
	((float)params_a(1) + (float)params_a(2) * Integral::cos(arg1) +
	 ((float)params_a(2) - (float)params_a(1)) * Integral::cos(arg2));
    }
  }
  
  // exit gracefully
  //
  return true;
}

// method: generateDolphChebyshev
//
// arguments: 
//  VectorFloat& output: (output) the window function
//  const VectorFloat& params: (input) window function parameters
//  long num_points: (input) number of points in the window
//
// return: a boolean value indicating status
//
// this method generates a Dolph-Chebyshev window function:
//
//  w(n) = { G/(2M+1) * [1/gamma + 
//         {             2 * sum Cheb(beta*k*PI/(2M+1)) *
//         {             cos(2*n*k*PI/(2M+1))]                  n <= |M|
//         { 0                                                  elsewhere
//
//  there are a couple of important auxiliary equations:
//
//   alpha_s = (2*M * 2.285 * (float)dw + 16.4) / 2.056
//
//  where dw is the transition bandwidth in normalized frequency.
//  using alpha_s, we compute gamma and beta (needed above):
//
//   gamma = 10 ** (-alpha_s / 20.0)
//   beta = cosh((1/2M) * acosh(1/gamma))
//
// a good reference for this window can be found at:
//
//  S.K. Mitra, Digital Signal Processing,
//  McGraw-Hill, Boston, Massachusetts, USA, 2001, pp. 456.
//
// there are two parameters allowed for this window:
//
//  params(0): G, the gain of the window
//  params(1): dw, the normalized transition bandwidth
//
// for an odd number of points, the window is defined as above.
// for an even number of points, an interpolation technique is used.
//
boolean Window::generateDolphChebyshev(VectorFloat& output_a,
				       const VectorFloat& params_a,
				       long num_points_a) const {
  
  // compute some constants related to the length of the window
  //
  long N = num_points_a;
  long Nm1 = N - 1;
  long Nm2 = Nm1 - 1;
  float scale = (float)Nm2 / (float)Nm1;
  
  // compute the half-window length
  //
  long M = (N - 1) / 2;
  
  // create output space
  //
  output_a.setLength(num_points_a);
  
  // in this approach, the constant represents the transition bandwidth.
  // the window length is supplied by the user. we compute the appropriate
  // constants to satisfy this relationship.
  //
  float A = (2*M * 2.285 * (float)params_a(1) + 16.4) / 2.056;
  float gamma = Integral::pow(10.0, -A / Integral::DB_MAG_SCALE_FACTOR);
  float gamma_inv = 1.0 / gamma;
  
  // from gamma and M, we compute beta
  //
  float beta = Integral::cosh(Integral::acosh(gamma_inv) / (float)(2*M));
  
  // case 1: an odd number of points
  //
  if ((N % 2) != 0) {
    for (long n = -M; n <= M; n++) {
      
      // compute the inner summation
      //
      float sum = 0;
      for (long k = 1; k < M; k++) {
	float arg1 = (float)k * Integral::PI / (float)(N);
	float arg2 = 2 * n * arg1;
	float val;
	Chebyshev::compute(val, beta * Integral::cos(arg1), k);
	sum += val * Integral::cos(arg2);
      }
      
      // save the output value
      //
      output_a(n + M) = params_a(0) / (float)(N) * (gamma_inv + 2 * sum);
    }
  }
  
  // case 2: an even number of points
  //  use a resampling strategy
  //
  else {
    for (long n = 0; n < N; n++) {
      float scaled_index = -M + scale * n;
      
      // compute the inner summation
      //
      float sum = 0;
      for (long k = 1; k < M; k++) {
	float arg1 = (float)k * Integral::PI / (float)(Nm1);
	float arg2 = 2 * scaled_index * arg1;
	float val;
	Chebyshev::compute(val, beta * Integral::cos(arg1), k);
	sum += val * Integral::cos(arg2);
      }
      
      // save the output value
      //
      output_a(n) = params_a(0) / (float)(Nm1) * (gamma_inv + 2 * sum);
    }
  }
  
  // exit gracefully
  //
  return true;
}

// method: generateGaussian
//
// arguments: 
//  VectorFloat& output: (output) the window function
//  const VectorFloat& params: (input) window function parameters
//  long num_points: (input) number of points in the window
//