자체 개발한 Digital Band Pass Filter 함수

개발이라고 할 것까진 없지만 간단하게 원하는 주파수 부분만 남기고 모두 제거하는 함수를 만들어 봤습니다.

 

arm_cfft_f32 함수로 시간 영역(Time Domain)인 샘플링 신호를 주파수 영역(Frequency Domain)으로 변환할 수 있었습니다. 여기서 10KHz 성분만 추출하여 다시 역 FFT를 수행하였습니다. 특정 주파수 영역만 추출하기 위해 아래와 같은 Band Pass Filter(BPF) 함수를 따로 만들었습니다. 참고로 다음 포스트에서 살펴보겠지만 CMSIS DSP 라이브러리는 FIR 필터 함수도 제공합니다. 이 필터 함수에 계수값만 조정하면, arm_cfft_f32 함수와 자체 제작한 BPF 함수를 사용할 때처럼 시간 영역과 주파수 영역을 오갈 필요 없이, 바로 시간 영역에서 시간 영역으로 원하는 결과만 필터링할 수 있습니다.

 

arm_cfft_f32 함수와 마찬가지로 사용자가 샘플링 주파수를 미리 알고 있다는 것을 전제로 하고 (fft_size-1) 인덱스가 샘플링 주파수를 의미하며 원하는 Pass 주파수 영역을 의미하는 bin 인덱스를 cutoff_freq_index 파라미터에 입력해야 합니다. band는 cutoff_freq_index 인덱스 양 옆으로 몇 인덱스까지 포함하는지를 나타냅니다.

 

/* @brief delete bins except for cutoff_freq_index ± band and mirror image
 *
 * @param cfft_result			result of complex fft including imaginary number
 * @param cutoff_freq_index		converted to cutoff frequency(= f_sampling * (cutoff_freq_index / (fft_size-1) )
 * @param band					width of cutoff frequency
 * @fft_size					number of bins(= Real Number + Imaginary Number)
 *
 * @return						arm_status
 *
 * */
static arm_status bandpass_filter(float32_t* cfft_result, uint32_t cutoff_freq_index, uint32_t band, uint32_t fft_size) {
	uint32_t nz1_max_index = (fft_size / 2) - 1;
	uint32_t nz1_start_index;
	uint32_t nz1_end_index;

	// mirror = NZ#2
	uint32_t mirror_index;
	uint32_t mirror_max_index = fft_size - 1;
	uint32_t mirror_start_index;
	uint32_t mirror_end_index;

	// check parameters
	if(cutoff_freq_index > (fft_size / 2)) { // overflowed NZ#1
		return ARM_MATH_ARGUMENT_ERROR;
	}

	mirror_index = fft_size - cutoff_freq_index;

	if(band > cutoff_freq_index) {
		nz1_start_index = 0; // bin min index
		mirror_end_index = mirror_max_index;
	} else {
		nz1_start_index = cutoff_freq_index - band;

		if((mirror_index + band) > mirror_max_index) {
			mirror_end_index = mirror_max_index;
		} else {
			mirror_end_index = mirror_index + band;
		}

	}

	if((cutoff_freq_index + band) > nz1_max_index) {
		nz1_end_index = nz1_max_index;
		mirror_start_index = nz1_max_index + 1; // mirror min index
	} else {
		nz1_end_index = cutoff_freq_index + band;

		if((mirror_index - band) < (nz1_max_index + 1)) {
			mirror_start_index = nz1_max_index + 1; // mirror min index
		} else {
			mirror_start_index = mirror_index - band;
		}
	}

	// NZ #1
	for(int i=0; i<(fft_size/2); i++) {
		if((i < nz1_start_index) || (i > nz1_end_index)) {
			cfft_result[i*2] = 0;
			cfft_result[i*2 + 1] = 0;
		}
	}

	// NZ #2
	for(int i=(fft_size/2); i<fft_size; i++) {
		if((i < mirror_start_index) || (i > mirror_end_index)) {
			cfft_result[i*2] = 0;
			cfft_result[i*2 + 1] = 0;
		}
	}

	return ARM_MATH_SUCCESS;
}

 

bandpass_filter 함수를 상용해 10KHz 주파수 성분을 제외한 나머지 주파수 성분을 모두 제거하고, 역FFT를 수행했습니다. 결과는 다음과 같습니다. 다음 포스트에서는 매트랩과 CMSIS FIR 함수를 사용해 BPF를 만들어 보고 제가 CMSIS FFT 함수와 연계하여 만든 BPF 함수와 결과와 속도를 비교해 보겠습니다.