formant.c

speech signal process tools

/* formant.c */
/*
 * This material contains unpublished, proprietary software of 
 * Entropic Research Laboratory, Inc. Any reproduction, distribution, 
 * or publication of this work must be authorized in writing by Entropic 
 * Research Laboratory, Inc., and must bear the notice: 
 *
 *    "Copyright (c) 1987-1990  AT&T, Inc.
 *    "Copyright (c) 1986-1990  Entropic Speech, Inc. 
 *    "Copyright (c) 1990-1994  Entropic Research Laboratory, Inc. 
 *                   All rights reserved"
 *
 * The copyright notice above does not evidence any actual or intended 
 * publication of this source code.     
 *
 * Written by:  David Talkin
 * Revised by: John Shore
 *
 */

static char *sccs_id = "@(#)formant.c	1.23	11/5/96	ATT/ESI/ERL";

#define VERSION "1.19"
#define DATE "9/20/94"

#include <Objects.h>
#ifndef LINUX
#include <sgtty.h>
#endif
#ifndef APOLLO_68K
#include <math.h>
#else
#include "math.h"
#endif
#include "tracks.h"
#include <esps/esps.h>
#include <esps/fea.h>
#include <esps/sd.h>
#include <string.h>

#define SYNTAX \
USAGE("formant [-p preemphasis] [-r n:[[+]m]] [-n num_formants] [-o lpc_order]\n \
 [-O output_path] [-i frame_interval] [-w window_duration] [-W window_type]\n \
 [-t lpc_type] [-S] [-f ds_freq] [-B max_buff_bytes] [-M maxrms_value]\n \
 [-N nominal_f1_freq] [-R maxrms_duration] [-x debug_level] [-F] infile.sd")

char *get_cmd_line();
void e_make_sd_hdr();
void e_make_fea_hdr();

int debug_level = 0, cmap_size;
int debug = 0;
int w_verbose = 0;

long max_buff_bytes = 2000000;  /* limit on size of signal data buffers */

int read_poles();

Signal *dpfund(), *lpc_poles(), *downsample(), *highpass(), *get_signal(),
       *dpform(), *new_signal(), *get_signal_from_path();

extern int fea_sd_special; /*in header.c*/
char *strcat(), *strrchr(), *strncpy(), *localloc();

char *new_ext(oldn,newex)
char *oldn, *newex;
{
  int	j;
  char *dot;
  static char newn[256];
  static int len = 0;

  dot = strrchr(oldn,'.');
  if(dot != NULL){
    *dot = 0;
    j = strlen(oldn) + 1;
    *dot = '.';
    
  }
  else j = strlen(oldn);
  if((len = j + strlen(newex)+ 1) > 256) {
    printf("Name too long in new_ext()\n");
    return(NULL);
  }
  strncpy(newn,oldn,j);
  newn[j] = 0;
  return(strdup(strcat(newn,newex)));
}
  
/*      ----------------------------------------------------------      */
main(argc,argv)
int	argc;
char	**argv;
{
  POLE **pole;
  int nframes, nform, nsamp, nsamp_ds, i,j, lpc_ord, in, lpc_type, f0only,
  w_type, save_sig_files, c;
  double  frame_int, wdur, maxrms_dur, maxrms_val, f0typ, 
  **freq, **band, ds_freq, nom_f1 = -10.0,
  preemp;
  double cor_wdur, fmin, fmax, *f0, **dpp, **dpps;
  short *data, *scr, *vuv;
  CROSS	**cp;
  DPREC	**dp;
  char *polename, *f0name, *sp, station[200], path[200];
  Signal *polesig, *f0sig, *formantsig, *spsig, *dssig, *hpsig, *stasig;
  Signal *hpsrc, *polesrc, *f0src;
  extern int optind;
  extern char *optarg;
  char *command;		/* function command line */
  char *range;
  struct header *ehd_sp, *ehd_ds, *ehd_hp, *ehd_fb, *ehd_f0;
  /* ESPS headers */

  /* make sure FEA_SD are read as shorts, since some routines insist on 
     them; note that this can be a problem if the external file does not 
     contain shorts since it's likely that we'll lose accuracy 
     */
	
  fea_sd_special = 1;

  command = get_cmd_line(argc, argv);
  lpc_ord = 12;
  lpc_type = 0;			/* use bsa's stabilized covariance if != 0 */
  w_type = 2;			/* window type: 0=rectangular; 1=Hamming; 2=cos**4 */
  ds_freq = 10000.0;
  wdur = .049;			/* for LPC analysis */
  cor_wdur = .01;		/* for crosscorrelation F0 estimator */
  maxrms_dur = 0.0;		/* window for max rms - 0 means all of past */
  maxrms_val = 0.0;		/* value for max rms - 0 means normalize */
  frame_int = .01;
  preemp = .7;
  nform = 4;
  fmin = 60.0;			/* limits for F0 estimator */
  fmax = 500.0;
  f0only = FALSE;
  save_sig_files = FALSE;
  *station = 0;			/* "stationarity" signal */
  *path = 0;
  range = NULL;
  if(argc < 2) 
    SYNTAX;

  while((c = getopt(argc,argv,"Ff:x:B:o:O:W:t:p:w:i:r:Ss:n:R:M:N:y:z:")) != EOF) {
    switch(c) {
    case 'F':
      f0only = TRUE;
      break;
    case 'f':
      ds_freq = atof(optarg);
      break;
    case 'x':
      debug_level = debug = atoi(optarg);
      break;
    case 'N':
      nom_f1 = atof(optarg);
      break;
    case 'o':
      lpc_ord = atoi(optarg);
      break;
    case 'O':
      strcpy(path,optarg);
      break;
    case 'W':
      w_type = atoi(optarg);
      break;
    case 't':
      lpc_type =  atoi(optarg);
      break;
    case 'p':
      preemp = atof(optarg);
      break;
    case 'w':
      wdur = atof(optarg);
      break;
    case 'i':
      frame_int = atof(optarg);
      break;
    case 'y':
      fmax = atof(optarg);
      break;
    case 'z':
      fmin = atof(optarg);
      break;
    case 'r':
      range = optarg;
      break;
    case 'S':
      save_sig_files = TRUE;
      break;
    case 's':
      strcpy(station,optarg);
      break;
    case 'n':
      nform = atoi(optarg);
      break;
    case 'B':
      max_buff_bytes = atol(optarg);
      break;
    case 'R':
      maxrms_dur = atof(optarg);
      break;
    case 'M':
      maxrms_val = atof(optarg);
      break;
    default:
      SYNTAX;
      break;
    }
  }

  /* Now process all of the input files
     (automatically creating output file names) */
  for(in=optind; in < argc; in++ ) {
    sp = argv[in];
    formantsig = polesig = f0sig = spsig = dssig = hpsig = NULL;

    /*
     * Check for errors in specifying parameters
     */
    if ((maxrms_val != 0.0) && (maxrms_dur != 0.0)) { 
      printf("Can't use both -R and -M\n");
      exit (1);
    }

    if ((maxrms_dur > 0.0) && (maxrms_dur < cor_wdur)) {
      maxrms_dur = frame_int;
      printf("WARNING: maxrms_duration less than cross-correlation window (%g)\n", cor_wdur);
      printf("\treset maxrms_duration to %g.\n", cor_wdur);
    }

    if(nform > (lpc_ord-4)/2){
      printf("Number of formants must be <= (lpc order - 4)/2 - exiting.\n");
      exit (1);
    }

    if(nform > MAXFORMANTS){
      printf("A maximum of %d formants are supported at this time - exiting.\n", 
	     MAXFORMANTS);
      exit(1);
    }
  
    if(range) {
      if((spsig = get_signal(sp,0.0,0.0,NULL))) {
	long start_p=0, end_p=0;
	double start, duration;
	lrange_switch (range, &start_p, &end_p, 1);
	if (end_p > spsig->file_size -1)
	    end_p = spsig->file_size -1;
	if(!end_p)
	  end_p = spsig->file_size -1;
	if(start_p > 0) start_p--; /* base 1 to base 0 conversion */
	if(start_p >= end_p) {
	  fprintf(stderr,"Unrealistic computation range limits\n");
	  exit(-1);
	}
	start = spsig->start_time + (((double)start_p)/spsig->freq);
	duration = ((double)(end_p - start_p + 1))/spsig->freq;
	spsig->utils->read_data(spsig,start,duration);
      }
    } else
      spsig = get_signal(sp,0.0,-1.0,NULL);

    if(spsig) {
      if(ds_freq != spsig->freq) {
	dssig = NULL;
	if(!((dssig = get_signal_from_path(path,new_ext(spsig->name,"ds"),0.0,-1.0,NULL)) &&
	     (dssig->freq == ds_freq) )) {
	  if(dssig) free_signal(dssig);
	  if((dssig = downsample(spsig,ds_freq))) {
	    if(spsig != dssig) {
	      /* insert proper ESPS header for downsampled signal*/
	      e_make_sd_hdr(dssig, spsig, command);
	      add_source_file(dssig->header->esps_hdr, sp, 
			      spsig->header->esps_hdr);
	      add_comment(dssig->header->esps_hdr, 
			  "This is the down sampled signal.\n");
	      dssig->header->magic = ESPS_MAGIC;
	      put_signal_to_path(path,dssig);
	    }
	  } else {
	    printf("Problems downsampling the signal %s\n",spsig->name);
	    exit(-1);
	  }
	}
      }

      hpsrc = (dssig ? dssig : spsig);

      if(f0only || (preemp < 1.0)) { /* be sure DC and rumble are gone! */
	if(! (hpsig = get_signal_from_path(path,new_ext(hpsrc->name,"hp"),0.0,-1.0,NULL))) {
	  if((hpsig = highpass(hpsrc))) {
	    /* insert proper ESPS header for high pass signal */
	    e_make_sd_hdr(hpsig, hpsrc, command);
	    add_source_file(hpsig->header->esps_hdr, hpsrc->name, 
			    hpsrc->header->esps_hdr);
	    add_comment(hpsig->header->esps_hdr, 
			"This is the high pass signal.\n");
	    hpsig->header->magic = ESPS_MAGIC;
	    put_signal_to_path(path,hpsig);
	  } else {
	    printf("Couldn't highpass the signal\n");
	    exit(-1);
	  }
	}
      }
    } else {
      printf("Can't access sampled speech signal: %s\n",sp);
      exit(-1);
    }
      
    polesrc = (hpsig ? hpsig : spsig);

    if(! f0only) {		/* compute LPC complex poles as needed  */
      if(!(polesig=lpc_poles(polesrc,wdur,frame_int,lpc_ord,preemp,lpc_type,
			     w_type))) {
	printf("Problems in lpc_poles()\n");
	exit(-1);
      }
      /* force writing of pole file as ASCII file*/
      polesig->header->magic = SIGNAL_MAGIC;
      /*	polesig->type |= SIG_ASCII;*/
      put_signal_to_path(path,polesig);
    }
      
    /*compute F0*/
      
    if(!(f0sig = dpfund(polesrc,cor_wdur, fmin, fmax, frame_int, 
			&cp, &dp, maxrms_dur, maxrms_val))) {
      printf("Problems in dpit()\n");
      exit(-1);
    }
    fflush(stdout);
    free_dp_cp(dp,cp,f0sig->buff_size); 
    /* insert proper header for F0 FEA file*/
    e_make_fea_hdr(f0sig, 5, command);
    f0sig->header->e_scrsd = 0;
    ehd_f0 = f0sig->header->esps_hdr;
    add_source_file(ehd_f0, polesrc->name, polesrc->header->esps_hdr);

    *add_genhd_d("preemphasis", (double *) NULL, 1, ehd_f0) = preemp;
    *add_genhd_d("window_duration", (double *) NULL, 1, ehd_f0) = wdur;
    *add_genhd_d("frame_duration", (double *) NULL, 1, ehd_f0) = frame_int;
    *add_genhd_s("lpc_order", (short *) NULL, 1, ehd_f0) = lpc_ord;
    *add_genhd_s("lpc_type", (short *) NULL, 1, ehd_f0) = lpc_type;
    *add_genhd_s("window_type", (short *) NULL, 1, ehd_f0) = w_type;
    (void) add_fea_fld("F0", 1L,
		       0, (double *) NULL, DOUBLE, (char **) NULL, ehd_f0);

    (void) add_fea_fld("prob_voice", 1L,
		       0, (double *) NULL, DOUBLE, (char **) NULL, ehd_f0);

    (void) add_fea_fld("rms", 1L,
		       0, (double *) NULL, DOUBLE, (char **) NULL, ehd_f0);

    (void) add_fea_fld("ac_peak", 1L,
		       0, (double *) NULL, DOUBLE, (char **) NULL, ehd_f0);

    (void) add_fea_fld("k1", 1L,
		       0, (double *) NULL, DOUBLE, (char **) NULL, ehd_f0);
    f0sig->header->magic = ESPS_MAGIC;
    put_signal_to_path(path,f0sig);
    f0sig->header->magic = SIGNAL_MAGIC;
    f0sig->type = P_DOUBLES | SIG_F0;