ft_dst_0.cc

这是处理语音信号的程序

// file: $PDSP/class/fourier_transform/v3.0/ft_dst_0.cc
//

// system include files
//
#include <math.h>

// isip include files
//
#include "fourier_transform.h"
#include "ft_qf.h"

#include "fourier_transform_constants.h"

// method: qf_real_dst_cc
//
// This code very closely follows the algorithm described in:
//     H.Guo et.al,"The Quick Discrete Fourier Transform", IEEE Transactions
//    on Acoustics Speech and Signal Processing 1994, pp. 445-448, 1994 
//
// This implementation is  adapted from the code written by G.A.Sitton at
// Rice University.
//
// arguments:
//  float_8* coeff_a : output dst coefficient array
//  float_8* data_a  : input data array
//
// return: a logical_1 value indicating status
//
// implements a disrete sine transform
//
logical_1 Fourier_transform::qf_real_dst_cc(float_8* coeff_a,
						     float_8* data_a) {
  // declare local variables
  //
  float_8  t0, t1, t2;
  float_8* XX;
  float_8* xx ;
  int_4	i, j, k, n, nq, n2, n2m1, temp;
  
  // Increment level & start	
  //
  n = nc[++ii];
  xx = ws[ii];
  
  // Last recursion level ?	
  //
  if (n == N2) {

    // Form odd function; eqn (24)
    //
    xx[0] = data_a[0] - data_a[6];
    xx[1] = data_a[1] - data_a[5];
    xx[2] = data_a[2] - data_a[4];
    
    // Save odd DST coefficients
    //
    t0 = (xx[0] + xx[2]) * sqrt2;
    
    coeff_a[1] = t0 + xx[1];
    coeff_a[3] = xx[0] - xx[2] ;
    coeff_a[5] = t0 - xx[1];
    
    // Form twiddled even function; eqn (24)
    //
    xx[0] = (float_8)(data_a[0] + data_a[6]) * sec1;
    xx[1] = (float_8)(data_a[1] + data_a[5]) * sec2;
    xx[2] = (float_8)(data_a[2] + data_a[4]) * sec3;
    
    // Save even DST coefficients
    //
    t0 = (xx[0] + xx[2]) * sqrt2;
    coeff_a[0] = (t1 = (t0 + xx[1])) + data_a[3];
    coeff_a[2] = (t2 = (xx[0] - xx[2])) + t1 - data_a[3];
    coeff_a[4] = (t0 = (t0 - xx[1])) + t2  + data_a[3];
    coeff_a[6] = t0 - data_a[3];
  }
  else {
    temp = n;
    n2 = temp >> 1;
    n2m1 = n2 - 1;
    nq = ic[ii];
    
    if ((k = mc[ii]) > n2m1)
      k = n2m1;
    
    // Form even-odd functions
    //
    if (nn <= n2) {
      temp = nn - 1;
    for (i = 0, j = nq; i < temp; ++i, j += nq)
      xx[i + n2m1] = (xx[i] = (float_8)data_a[i]) * sc[j];
    }
    else  {
      temp = n - nn;
      for (i = 0, j = nq; i < temp; ++i, j += nq)
	xx[i + n2m1] = (xx[i] = (float_8)data_a[i]) * sc[j];
      
      for (;i < n2m1; ++i, j += nq){
	xx[i] = data_a[i] - data_a[(n - 2) - i];
	xx[i + n2m1] = ((float_8)data_a[i] +
			(float_8)data_a[(n - 2) - i]) * sc[j];
      }
    }
    
    // Do recursive DSTs; similar to eqn (26)
    //
    qf_real_dst_cc((float_8*)(XX = &xx[n - 2]), (float_8*)xx);
    qf_real_dst_cc((float_8*)&XX[n2m1], (float_8*)&xx[n2m1]);
    
    // Save DST coefficients
    //
    coeff_a[1] = XX[0];
    
    if (n2 >= nn){
      t0 = coeff_a[0] = XX[n2m1];
      
      for (i = 1; i < k; i += 2) {
	coeff_a[j = i + i] = (t1 = XX[i + n2m1]) + t0;
	coeff_a[j + 1] = XX[i];
	coeff_a[j + 2] = (t0 = XX[i + n2m1 + 1]) + t1;
	coeff_a[j + 3] = XX[i + 1];
      }
      coeff_a[i + i] = t0;
    }
    else {
      coeff_a[0] = (t1 = XX[n2m1]) + (t0 = data_a[n2m1]);
      
      for (i = 1; i < k; i += 2){
	coeff_a[j = (i + i)] = (t2 = XX[i + n2m1]) + t1 - t0;
	coeff_a[j + 1] = XX[i];
	coeff_a[j + 2] = (t1 = XX[i + n2m1 + 1]) + t2 + t0;
	coeff_a[j + 3] = XX[i + 1];
      }
      
      coeff_a[i + i] = t1 - t0;
    }
  }
  
  // Decrement level & exit
  //
  --ii;
  
  // exit gracefully
  //
  return ISIP_TRUE;
  
}