看板 Math 關於我們 聯絡資訊
直接寫程式來湊答案 列出主程式 #define M_PI 3.14159265358979323846 #define SAMPLE_RATE 4 #define FFT_N 4 complex<double> s[FFT_N] = {0}; void test() { int i; int sampleRate = SAMPLE_RATE; //Hz double deltaDeg = 2 * M_PI / sampleRate; TRACE("input:\n"); for(i = 0; i < sampleRate; i++) { double x = deltaDeg * i; s[i] = 1 + sin(x);// + 2 * sin(2 * x) + 3 * sin(3 * x); TRACE("s[%d]=%f\n", i, s[i]); } FFT(s, FFT_N); TRACE("\nafter FFT\n"); for(i = 0; i < FFT_N; i++) { double real = s[i].real(); double imag = s[i].imag(); TRACE("s[%d]=(%f,%f) = %f deg %f\n", i, real, imag, sqrt(real * real + imag * imag), atan(imag / real) / M_PI * 180 ); } } -- ※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 111.82.254.91 (臺灣) ※ 文章網址: https://www.ptt.cc/bbs/Math/M.1571655907.A.7B7.html 以上程式 s[i] = 1 + sin(x); 因為用的是 Visual C++ 6.0 with MFC 所以 s[i] 這種寫法,是直接填入實部 我用的測試函數是含直流成份 1, 1HZ sin 波,相角 0 (超級單純) 在此我們以 4HZ 做取樣,填入 size 為 4 的陣列中 #define SAMPLE_RATE 4 #define FFT_N 4 來看輸出 input: s[0]=1.000000 s[1]=2.000000 s[2]=1.000000 s[3]=0.000000 after FFT s[0]=(4.000000,0.000000) = 4.000000 deg 0.000000 s[1]=(-0.000000,-2.000000) = 2.000000 deg 90.000000 s[2]=(0.000000,0.000000) = 0.000000 deg -1.#IND00 s[3]=(-0.000000,2.000000) = 2.000000 deg -90.000000 s[0] 是直流成份,只有實部,其值為 4 s[1] 是 1hz, 只有虛部,其值為 -2 這什麼意思? https://www.itread01.com/content/1544320336.html 》在輸出訊號而言,相鄰兩個取樣點的間隔為頻率解析度fs/N, 》其中fs為取樣頻率,N為輸入訊號的取樣點數 所以頻率解析度 = 4 / 4 = 1 果然如此,是這樣解讀頻域的 XD 但有沒有文件能告訴我,為什麼 s[0] 要除以 4,才變成直流成份 1 而 s[1] 要除以 2,才變成振幅 1 而且為什麼它是表達成 (0, -2),這是 -90度喔! 但我明明原函式是 0度 (程式運算 atan 時,在 90 附近會有誤差,差個負號,結果角度錯了) ※ 編輯: HuangJC (111.82.254.91 臺灣), 10/21/2019 19:28:48 我們改用 s[i] = 1 + sin(x + 30.0 / 180 * M_PI); 相角 30 度當成輸入波型去做好了 after FFT s[0]=(4.000000,0.000000) = 4.000000 deg 0.000000 s[1]=(1.000000,-1.732051) = 2.000000 deg -60.000000 s[2]=(0.000000,0.000000) = 0.000000 deg -1.#IND00 s[3]=(1.000000,1.732051) = 2.000000 deg 60.000000 就可以看到 s[1] 輸出是 -60 度,也就是再加 90度就可以還原原波型 為什麼要這樣迂迴啊... ※ 編輯: HuangJC (111.82.254.91 臺灣), 10/21/2019 19:31:04 好,那我現在來驗算 取樣 500 點,因為必需補到 2^N 所以陣列大小是 512 個,可以修改我的程式變成這樣來測 #define SAMPLE_RATE 500 #define FFT_N 512 還是用單純的波型 s[i] = 1 + sin(x); after FFT s[0]=(500.000000,0.000000) = 500.000000 deg 0.000000 s[1]=(6.641283,-253.006292) = 253.093442 deg -88.496361 s[2]=(-13.022571,6.307806) = 14.469823 deg -25.844379 s[3]=(-12.571770,1.665885) = 12.681663 deg -7.548286 輸出結果,直流成份要 除以 500 那麼 s[1] 的值除以 250,大約是 1 不過我們計算到 253,角度也不再是 -90 度 換言之 500 補到 512 還是有影響的 這東西只是'可參考',但補 0 不是不影響原波型的! ※ 編輯: HuangJC (111.82.254.91 臺灣), 10/21/2019 19:44:18
HuangJC : 接下來更挑戰了 Orz,明天再找同學確認 10/22 01:31
HuangJC : 這套網路上抓來的 FFT,基本上處理最多 512點 10/22 01:32
HuangJC : 可是我收到的檔案快六萬點,這怎麼處理! 10/22 01:32