DOCUMENTATION FOR FFT.C Documentation dated 10/89. Program and documentation Copyright (c) 1989, E.R.I. All rights reserved. INTRODUCTION FFT.C is the source file for supporting common FFT operations. All FFT operations in FFT.C use double-precision floating-point data. FFT.C was compiled under Microsoft C, v. 5.1, using the emulation library. Thus, it will use a math coprocessor if one is installed on the computer, but does not require it. As with all my compilations, the /Od flag was set, disabling optimization. FFT.C includes one function that reads data files (data_multiply()). The data format used in data_multiply() files is: A B C data data data; where A is the number of data points in the file, B is the minimum value of the data set, and C is the maximum value in the data set. Each data point must be separated by white space, which can but does not have to include a carriage return. The data files must be text format, not binary. FFT.C makes explicit calls to functions in MENU.C and MOUSELIB.C. MENU.C is part of FTPLOT.LIB, and MOUSELIB.C forms a separate library, MOUSELIB.LIB. EXTERNALLY-VISIBLE FUNCTIONS IN FFT.C The FFT.C functions are declared in FFT.H, which should be included in any program that utilizes these routines. FUNCTION: data_multiply() TYPE: void ARGUMENTS: Two character strings data_multiply(file1, file2) multiplies the data in file1 by the data in file2, and outputs the result to a file that the user specifies in response to a prompt. The format for each of these files is that explained above. FUNCTION: dbl_scale() TYPE: int ARGUMENTS: three doubles and two ints i=dbl_scale(val,min,max,out_min,out_max) scales a double-precision number VAL, that can range from MIN to MAX into an int. The output int can range from OUT_MIN to OUT_MAX. dbl_scale() is used to scale a double-precision number for graphing. Notice that you can make OUT_MIN larger than OUT_MAX if you want to invert the scale for screen display -- that is, most computer screens have large Y values at the BOTTOM and small values at the top. By making OUT_MIN the largest screen y-value and OUT_MAX the smallest, the routine will auto- scale with inversion. FUNCTION: fft() TYPE: void ARGUMENTS: two arrays of doubles and two ints. fft(x,y,n,flag) accepts the n data points x[i], y[i] and transforms them to the other domain. n must be an even power of 2, and cannot exceed 4096. If flag is -1, then x and y are time-domain data points (x=real, y is normally an array of 0s). If flag is 1, then x and y are frequency-domain data points (rectangular coordinates, x=real, y=imaginary). The routine is an in-place algorithm, which means that the result of the fft is returned in the x and y arrays. The fft routine does not normalize forward (time to frequency) transforms. Rather, it waits until the inverse transform to normalize. If you perform a forward and then an inverse transform, you will end up with the original data (within round-off error). However, if you only perform forward transforms, the data will be weighted by the number of data points in the transform. For example, if you do a 256-point fft of a 1 Hz sine wave, and a 512-point fft of the same 1 Hz sine wave, the unnormalized data will make the amplitude of the 512-point transform appear to be twice as large as the 256- point transform. If you are only comparing transforms with the same number of data points, the failure to normalize does not matter in either domain. However, if you will be comparing different size transforms, you should normalize the results by dividing each transformed data point by the number of points in the data set. FUNCTION: get_filter() TYPE: pointer to double ARGUMENTS: one array of doubles, two ints, and four doubles. fptr = get_filter(filter, n, direction, nyquist, cutoff, rolloff, units) builds a high-pass or low-pass filter, and returns a pointer to the array of frequency-domain filter coefficients (fptr). filter is the array that the coefficients will be written into, so fptr = &filter[0]. n is the number of points in the filter. direction = 1 specifies that the filter is high-pass, and direction = -1 specifies that the filter is low-pass. nyquist is the maximum defined frequency in the filter data set. cutoff is the cutoff frequency of the filter, rolloff is the rolloff in either dB per octave or dB per decade, depending on the value of units (2.0 or 10.0). get_filter() only produces positive/negative format amplitude data. If you want positive-only data, throw away the last n/2 data points. The filter generated by get_filter() is used by multiplying an amplitude spectrum by filter. If you save filter as a file in the format given in this document's introduction, you can use data_multiply() on the amplitude part of polar-format frequency data to perform the filtering. FUNCTION: int_scale() TYPE: int ARGUMENTS: five ints. i=int_scale(val,min,max,out_min,out_max) is the integer-input version of dbl_scale(). See that function for details. FUNCTION: long_scale() TYPE: int ARGUMENTS: three long ints and two ints. i=long_scale(val,min,max,out_min,out_max) is the long int version of dbl_scale(). See that function for details. FUNCTION: minmax() TYPE: void ARGUMENTS: one array of doubles, one int, and two pointers to doubles minmax(array, n, minv, maxv) finds the minimum and maximum value in the n data points in array, and writes those values into the locations specified by minv and maxv. This is useful in generating the data for the file format specified in the introduction. FUNCTION: polar_to_xy() TYPE: void ARGUMENTS: four arrays of doubles and one int. polar_to_xy(amp, phase, x, y, n) converts the n polar-format data points in amp, phase to rectangular-format data in x, y. amp,phase and x,y can be the same arrays if you want to overwrite the current data, i,e,, polar_to_xy(x,y,x,y,n) a valid call. The data arrays in fft() are assumed to be x, y format. FUNCTION: unscramble_transform() TYPE: void ARGUMENTS: one array of doubles and one int. unscramble_transform(array, n) takes the positive/negative format frequency data array array and places the negative frequency components at the beginning of the array, so the x-axis of the graph is continuous. When you perform a forward fft, the data array is ordered as follows: freq 0, freq 1, freq 2, ...freq nyquist, [the negative freq one less than the nyquist freq, the negative freq 2 less than the nyquist freq], etc. I unscramble the array by moving freq 0 up to BUT NOT INCLUDING freq nyquist to the end of the array. The x-axis thus runs from -freq nyquist to +freq (nyquist-1), in order. The reason that I don't move the Nyquist frequency is that two successive calls to this function will return the array to its original order (which you must do before performing an inverse transform on the data set). FUNCTION: xy_to_polar() TYPE: void ARGUMENTS: four arrays of doubles and an int. xy_to_polar(x, y, amp, phase, n) converts the n data points in rectangular- format frequency arrays x and y to polar-format data arrays amp and phase. If you do not want to preserve the original format data, you can use the same data names for x and amp and y and phase, i.e., xy_to_polar(x, y, x, y, n) is a legal call.