filter.c

来自「speech signal process tools」· C语言 代码 · 共 816 行 · 第 1/2 页

/*
 * 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) 1986-1990  Entropic Speech, Inc. 
 *    "Copyright (c) 1990-1991  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:   Brian Sublett
 * Checked by:   Dave Burton
 * Revised by:   Bill Byrne (converted to FEA_FILT)
 *
 * This program reads data from an ESPS sampled-data file and
 * filters it with a filter of arbitrary length.  The resulting
 * output data is then written to another ESPS data file.  The
 * program may perform interpolation filtering, where the output
 * sampling rate is different from the input sampling rate.  The
 * program calls the routine block_filter to implement regular
 * filtering.  It uses interp_filt to implement interpolation
 * filtering.
 */

static char *sccs_id = "%W% %G% ERL"

#define Version "3.22"
#define Date "5/22/92"

THIS FILE NEEDS TO EDITED BEFORE 5.2 RELEASE TO FIX SKIPREC.  THIS STATEMENT
IS IN HERE TO CAUSE A COMPILE ERROR SO THAT WE DON'T FORGET THAT.

/* 
 * System Includes
 */
# include <stdio.h>
# include <math.h>

/*
 * ESPS Includes
 */
# include <esps/esps.h>
# include <esps/fea.h>
# include <esps/feasd.h>
# include <esps/feafilt.h>
# include <esps/unix.h>
# include <esps/limits.h>

# define AR_LIM 2048	/* Store this many outputs before
			    writing to a file.  */

# define FILNAMSIZE 40	/* Maximum filter name size.  */

# define SYNTAX \
USAGE("filter [-d data_type] [-f filter name] [-i up/down] [-{pr} range]\n [-x debug_level] [-z] [-F filt_file] [-P parfile] in_file out_file");

/*
 * ESPS Functions
 */

void lrange_switch();
char *get_cmd_line();
void put_sd_recd();
void add_genzfunc();
char *arr_alloc();


int debug_level = 0;

main(argc, argv)
    int	    argc;
    char    *argv[];
{
    FILE  *fpin = stdin, *fpout = stdout; /* in and out file stream ptrs*/
    FILE *fpco;				  /* strream ptr for FILT file*/
    struct header *ih, *oh, *fh = NULL;	 /* file header ptrs*/
    char *in_file = NULL, *out_file = "<stdout>"; /* In and Out file names*/
    int i, j, nsiz, dsiz, nn, up = 1, down = 1, outflag = 0, coutflag = 0, nout, cplxflag=0;
    int co_num = 1;
    struct feafilt *filt_rec;		/* FEAFILT file data structure */
    double  *x, *y, **state;     	/* state vector info */
    double  *cx, *cy, **cstate;     	/* state vector info */
    long    dim[2];			/* dimensions of state, cstate */
    float *fzeros, *fpoles;		/* to write zfunc */
    double *zeros, *poles;		/* holds zero and pole coeff */
    char co_source = 'p';		/* flag for source of coeffs */

    char *filt_file;			/* holds FILT filename */
    char na_name[FILNAMSIZE], da_name[FILNAMSIZE];  /* holds coeffs */
    char nsiz_name[FILNAMSIZE], dsiz_name[FILNAMSIZE]; /* number of coeff */
    short   type, interp = NO, fflag = NO; 
    int	    order;
    long    c, start = 1, nan = LONG_MAX, start_p = 0, end = LONG_MAX, n;

    struct zfunc *pzfunc, tmp_zfunc;
    char *param_file = NULL, *filter_name = "filter";

    struct feasd *irec;                 /* input record */
    struct feasd *orec;                 /* output record */
    double **dptr;                      /* input data pointer */
    double_cplx **cdptr;                /* input data pointer - complex */
    double **optr;                      /* output data pointer */
    double_cplx **coptr;                /* output data pointer - complex */

    int ch_no;                          /* channel number */
    int total_pts = 0;			/* number of output points*/    
    int num_of_files = 0;		/* number of input files*/
    int rflag = 0;			/* flag for command line range*/
    char *range;			/* holds range arguments*/
    int check_common;			/*flag used by read_params*/
    int dflag = 0;			/*flag for data type*/
    char data_type;			/*holds output data type*/
    long big_count = 0;			/* holds # values > maxvalue*/
    double maxvalue = 0;		/* set if -d option causes 
					   conversion to lower type*/
    int Fflag = 0;			/* indicates filt file input*/
    int zflag = 0;                      /* flag to do filter delay adj */
    double fabs();
    long num_channels;                  /* holds # of sampled data channels*/
    double src_sf;                      /* source and output sampling freq.*/
    extern char *optarg;
    extern optind;


/* Check the command line options. */

    while ((c = getopt(argc, argv, "P:f:F:x:d:i:p:r:z")) != EOF)
    {
	switch (c)
	{
	case 'P':
	    param_file = optarg;
	    break;
	case 'f':
	    filter_name = optarg;
	    fflag = YES;
	    break;
	case 'd':
	    dflag++;
	    data_type = optarg[0];
	    break;
	case 'F':
	    Fflag++;
	    filt_file = optarg;
	    co_source = 'f';
	    TRYOPEN("filter", filt_file, "r", fpco);
	    break;
	case 'x':
	    if (sscanf(optarg, "%d", &debug_level) != 1)
	    {
		Fprintf (stderr, "Error reading -x arg. Exiting.\n");
		exit(1);
	    }
	    break;
	case 'i':
	    if ((nn = sscanf(optarg, "%d/%d", &up, &down)) != 2)
	    {
		Fprintf(stderr, "filter:interpolation format is up/down\n");
		Fprintf(stderr, "filter: nn=%d up=%d down=%d\n", nn, up, down);
		exit(1);
	    }
	    if (up <= 0)
	    {
	    Fprintf(stderr, "filter: up = %d is illegal; must be > 0.\n", up);
	    exit(1);
	    }
	    if (down <= 0)
	    {
	    Fprintf(stderr, "filter: down = %d is illegal; must be > 0.\n", down);
	    exit(1);
	    }
	    interp = YES;
	    break;
	case 'p':
	case 'r':
	    rflag = 1;
	    range = optarg;
	    break;
        case 'z':
	    zflag++;
	    break;
	default:
	    SYNTAX;
	}
    }

/* Change the filtername default if -F was used and -f wasn't. */

    if (fflag == NO && co_source == 'f')
	filter_name = "1";

/* Get the filenames. */

/* 
 * Get the number of filenames on the command line
 */ 
    if ((num_of_files = argc - optind) > 2){
	Fprintf(stderr, "Only two file names allowed\n");
	SYNTAX;
    }
    if (debug_level > 0)
	Fprintf(stderr, "filter: num_of_files = %d\n", num_of_files);

    if (num_of_files == 0){
	Fprintf(stderr, "No output filename specified.\n");
	SYNTAX;
    }


    if (num_of_files == 2) {

/*
 * Both input and output file names specified on command line
 */
	in_file = eopen("filter", argv[optind++], "r", FT_FEA, FEA_SD,
			&ih, &fpin);
	if (debug_level) 
	    Fprintf(stderr, "filter: Input file is %s\n", in_file);


    	out_file = eopen("filter", argv[optind], "w", NONE, NONE,
			    (struct header **)NULL, &fpout);
	if (debug_level) 
	    Fprintf(stderr, "filter: Output file is %s\n", out_file);

	/*
	 * Decide whether to check COMMON or not
	 */
	if (strcmp(in_file, "<stdin>") == 0)
	    check_common = SC_NOCOMMON;
	else
	    check_common = SC_CHECK_FILE;
	
    }

    if ( num_of_files == 1 )
        check_common = SC_CHECK_FILE;

    /* Read the parameter file. */
    (void)read_params(param_file, check_common, in_file);

    if (num_of_files == 1)  {
/*
 * Only output file specified on command line
 */
	if (in_file == NULL){
	    if (symtype("filename") == ST_UNDEF){
	        Fprintf(stderr, 
			"filter: No input file on command line or in Common\n");
		exit(1);
	    }
		
	    in_file = eopen("filter", getsym_s("filename"), "r",
			    FT_FEA, FEA_SD, &ih, &fpin);
	}
	
	if (strcmp( in_file, argv[optind] ) == 0 ) {
	    Fprintf(stderr, "filter: Fatal Error. Identical input and output names\n");
	    exit(1);
	}

	out_file = eopen("filter", argv[optind], "w", NONE, NONE,
			    (struct header **)NULL, &fpout);
	if (debug_level) 
	    Fprintf(stderr, "Output file is %s\n", out_file);
    }
/*
 * Check input types
 */
    if (is_field_complex(ih, "samples")){
	if (debug_level)
	    Fprintf(stderr, "filter: Complex data.\n");
	cplxflag = 1;
    }

    num_channels = get_fea_siz("samples", ih, (short *) NULL, (long **) NULL);
    if (debug_level) 
	Fprintf(stderr, "field samples is dimension %d\n", num_channels);

/* Create the output header. */

    oh = new_header(FT_FEA);

    add_source_file(oh, in_file, ih);
    add_comment(oh, get_cmd_line(argc, argv));

    Strcpy(oh->common.prog, "filter");
    Strcpy(oh->common.vers, Version);
    Strcpy(oh->common.progdate, Date);

  /*
   * Find output data type and set maximum value
   */
    if (dflag){
	switch (data_type)
	{
	case 'B':
	case 'b': 
	    maxvalue = CHAR_MAX;
	    if (cplxflag) 
		type = BYTE_CPLX;
	    else
		type = BYTE;
	    break;
	case 'S':
	case 's': 
	    maxvalue = SHRT_MAX;
	    if (cplxflag)
		type = SHORT_CPLX;
	    else
		type = SHORT;
	    break;
	case 'L':
	case 'l': 
	    maxvalue = LONG_MAX;
	    if (cplxflag)
		type = LONG_CPLX;
	    else
		type = LONG;
	    break;
	case 'F':
	case 'f': 
	    maxvalue = FLT_MAX;
	    if (cplxflag)
		type = FLOAT_CPLX;
	    else
		type = FLOAT;
	    break;
	case 'D':
	case 'd': 
	    maxvalue = DBL_MAX; 
	    if (cplxflag)
		type = DOUBLE_CPLX;
	    else
		type = DOUBLE;
	    break;
	default: 
	    Fprintf(stderr, 
		    "filter: -d type \"%c\" unknown - exiting.\n",
		    data_type);
	    exit(1);
	}
    }
    else{
        type = get_sd_type(ih);
	switch (type)
	{
	case CHAR: 
	case BYTE:
	case BYTE_CPLX:
	    maxvalue = CHAR_MAX;
	    break;
	case SHORT: 
	case SHORT_CPLX: 
	    maxvalue = SHRT_MAX;
	    break;
	case LONG: 
	case LONG_CPLX: 
	    maxvalue = LONG_MAX;
	    break;
	case FLOAT: 
	case FLOAT_CPLX: 
	    maxvalue = FLT_MAX;
	    break;
	case DOUBLE: 
	case DOUBLE_CPLX: 
	    maxvalue = DBL_MAX; 
	    break;
	default: 
	    Fprintf(stderr, "filter: input data type unknown\n");
	    exit(1);
	}
    }

    (void) add_genhd_d("max_value", &maxvalue, 1, oh);

/*
 * Initialize output fea_sd header
 */
    {
	double start_time = 0.;    /*this value is overwritten later*/

	src_sf = get_genhd_val("record_freq", ih, (double)0.);
	if (src_sf == 0){
	    Fprintf(stderr,
	        "filter: Input file has invalid \"record_freq\" - exiting.\n");
	    exit(1);
	}

        if (init_feasd_hd(oh, type, (int) num_channels, &start_time, 0, src_sf) 
	    != 0)
	{
	    Fprintf(stderr,
		"filter: Couldn't allocate fea_sd header - exiting\n");
	    exit(1);
	}

    }
       
/*
 * Get range 
 */

    if (!rflag){
	if (symtype("start") == ST_UNDEF){