ft_qf_2.cc

这是处理语音信号的程序

// file: $PDSP/class/fourier_transform/v3.0/ft_qf_2.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_init_cc
//
// This code very closely follows the algorithm described in:
//    H.Guo et.al,"The Quick Discrete Fourier Transform",  Proceedings of
//    ICASSP 1994, pp. 453-456, SanFrancisco, CA, April 1994.
//
// This implementation is  adapted from the code written by G.A.Sitton at
// Rice University.
//
// arguments:
//   int_4 order_a: (input) new order
//
// return: a logical_1 value indicating status
//
// creates lookup tables for sine and cosine terms. note that space is
// allocated in this method only when the order changes. also all temporary
// space is allocated here.

logical_1 Fourier_transform::qf_init_cc(int_4 order_a) {

  // if the order is the same as the previous call, we are done
  //
  if (qf_last_order_d == order_a) {
    return ISIP_TRUE;
  }
 
  // assign the order
  //
  N_d = order_a;
  qf_last_order_d = N_d; 

  float_8 q;
  int_4 i, m, n2, n2m1, temp;

  n2 = N_d >> (int_4)1;
  is_power_cc((int_4&)m, (int_4)2);
  
  // Allocate the QFT memory	
  //
  ws = new float_8*[sizeof(float_8*)*m + sizeof(int_4)*6*m +
		   sizeof(float_8)*3*(N_d + m)];

  // seperate memory is allocated for real and complex computations
  // in this case since the qft real routine is called from the complex
  // routine
  //
  qf_real_coeff_d = new float_8[N_d];
  qf_imag_coeff_d = new float_8[N_d];;
  qf_real_temp_d = new float_8[N_d];

  // memory space specifically used only for the complex routine
  //
  qf_comp_real_coeff_d = new float_8[2*N_d];
  qf_comp_imag_coeff_d = new float_8[2*N_d];;
  qf_comp_real_temp_d = new float_8[N_d];
  qf_comp_imag_temp_d = new float_8[N_d];
    
  // the following variables are used to do the indexing, recursion and
  // workspace traversing
  //
  nc = (int_4*)&ws[m];
  ic = &nc[m];
  mc = &ic[m];
  sc = (float_8*)&mc[m];
  ws[0] = &sc[n2 - 1];
  
  // Initialize tables	
  //
  nc[0] = N_d;
  ic[0] = 1;
  mc[0] = m;
  
  for (i = 1; i < m; ++i) {
    temp = nc[i - 1];
    nc[i] = temp >>  (int_4)1;
    temp = ic[i - 1];
    ic[i] = temp << (int_4)1;
    
    if (i > 1)
      ws[i] = ws[i - 1] + (2 * (int_4)(nc[i - 1] + 2));
    else
      ws[i] = ws[i - 1] + nc[i] + 1;	
  }	
  
  // Compute scaled secants  to be used to compute the even and odd functions
  // which are used in the DCT and DST computations.
  //
  q = M_PI/N_d;

  n2m1 = n2 - 1;
  for (i = 0; i < n2m1; ++i) {
    sc[i] = 0.5 / (float_8)cos(q * i);	
  }
  
  // exit gracefully
  //
  return ISIP_TRUE;
}