et-go/et_calc/algorithm/fft.go

package algorithm

import (
	"log"
	"math"
	"math/cmplx"
	"scientificgo.org/fft"
)

func IsPow2(N int) bool {
	if N == 0 {
		return false
	}
	return (uint64(N) & uint64(N-1)) == 0
}

// 频谱计算(全程分析)
func FFTforward(raw []float64, fs float64) (float64, []float64) {
	rawLen := len(raw)
	power := math.Ceil(math.Log(float64(rawLen)) / math.Log(2))

	signLen := int(math.Pow(2, power))

	signArray := make([]float64, signLen)
	copy(signArray[:], raw[:])

	if !IsPow2(signLen) {
		log.Panicf("fft源数据长度[%d]异常,不是2^n", signLen)
	}
	cpx := Float64ToComplex128Array(signArray)
	fftData := fft.Fft(cpx, false)
	fftLen := signLen / 2
	df := fs / float64(signLen)

	resp := make([]float64, fftLen)
	for i := 0; i < fftLen; i++ {
		resp[i] = cmplx.Abs(fftData[i]) * 2 / float64(signLen)
	}
	return df, resp
}

func dfftCalc(raw []float64) []float64 {
	cpx := Float64ToComplex128Array(raw)
	rp := fft.Fft(cpx, false)
	frp := Complex128ToFloat64Array(rp)
	return frp
}

func Float64ToComplex128Array(x []float64) []complex128 {
	y := make([]complex128, len(x))
	for i, v := range x {
		y[i] = complex(v, 0)
	}
	return y
}
func Complex128ToFloat64Array(x []complex128) []float64 {
	y := make([]float64, len(x))
	for i, v := range x {
		y[i] = real(v)
	}
	return y
}