// บทที่ 14.13 · อนุกรมฟูเรียร์ · From Fourier Series to DFT & FFT

จากอนุกรมฟูเรียร์สู่ DFT และ FFT: เครื่องวัดคำนวณฮาร์มอนิกอย่างไร

เครื่องวิเคราะห์กำลังไฟฟ้า ฟังก์ชัน FFT บนสโคป และเครื่องวัดพลังงานที่ทำเองด้วยไมโครคอนโทรลเลอร์ ไม่ได้อินทิเกรตแบบในตำรา มันสุ่มตัวอย่างแล้วรวมผลคูณ ซึ่งก็คืออินทิกรัลของสูตรออยเลอร์ในรูปผลบวก หน้านี้บอกว่าผลบวกนั้นตรงกับอนุกรมฟูเรียร์เมื่อไร และผิดได้อย่างไรบ้าง เพื่อให้อ่านตัวเลขบนเครื่องวัดได้อย่างเข้าใจ

เขียนโดย อาจารย์มดอัปเดต 12 ต.ค. 2569

อนุกรม → DFT อเลียสซิง Leakage IEC 61000-4-7 Goertzel + โค้ด

สรุปใน 30 วินาที

สุ่ม N จุดพอดีจำนวนเต็มคาบ: $c_k = X[k]/N$ · $A_k = 2|X[k]|/N$ · RMS $= \sqrt2|X[k]|/N$
ความละเอียด $\Delta f = f_s/N = 1/T_w$ · ฮาร์มอนิกที่เกิน $f_s/2$ พับลงมาทับฮาร์มอนิกต่ำ (อเลียสซิง)
หน้าต่างไม่ครบคาบ → leakage แก้ด้วยซิงก์ (PLL) หน้าต่าง Hann/flat-top หรือรวมกลุ่มช่อง
IEC 61000-4-7: 10 คาบ (200 ms) ช่องละ 5 Hz ฮาร์มอนิก h อยู่ช่อง 10h

ระดับพื้นฐาน ม.ปลาย · ปวส. · ปี 1กลาง ป.ตรี วิศวกรรมขั้นสูง ป.โท · งานวิจัย

จากอินทิกรัลสู่ผลบวก: DFT คืออนุกรมฟูเรียร์ของข้อมูลที่สุ่มมา กลาง

แทนอินทิกรัลใน C2 ด้วยผลรวมสี่เหลี่ยมผืนผ้าของ N ตัวอย่างที่ห่างกัน $T/N$ ได้สูตร DFT ทันที และถ้าสัญญาณมีฮาร์มอนิกไม่เกิน N/2 ผลรวมนี้ไม่ใช่ค่าประมาณ มันตรงทุกหลัก เพราะไซน์ที่สุ่มครบคาบยังตั้งฉากกันแบบไม่ต่อเนื่อง

D1 Discrete Fourier transform นิยาม DFT

X[k] = \sum_{n=0}^{N-1} x[n]\,e^{-j2\pi kn/N},\qquad k = 0, 1, \dots, N-1

เทียบกับ $c_k = \frac1T\int_T x(t)e^{-jk\omega_0 t}dt$ : แทน $t = nT/N$ , $dt \to T/N$ ได้ $c_k \approx X[k]/N$

D2 From DFT bins to harmonic values แปลงผล DFT เป็นค่าฮาร์มอนิก (หน้าต่างครบ m คาบพอดี)

\text{DC} = \frac{X[0]}{N},\qquad A_h = \frac{2|X[mh]|}{N},\qquad V_{h,rms} = \frac{\sqrt2\,|X[mh]|}{N},\qquad \varphi_h = \angle X[mh]

m = จำนวนคาบในหน้าต่าง ฮาร์มอนิกที่ h อยู่ช่อง $k = mh$ · ตรงทุกหลักถ้าไม่มีความถี่เกิน $f_s/2$

D3 Frequency resolution ความละเอียดความถี่

\Delta f = \frac{f_s}{N} = \frac{1}{T_w}

ขึ้นกับความยาวหน้าต่าง $T_w$ เท่านั้น ไม่ขึ้นกับอัตราสุ่ม อยากแยกความถี่ใกล้กันต้องเก็บข้อมูลนานขึ้น

D4 FFT คำนวณ DFT เร็ว

N^2 \;\to\; \tfrac{N}{2}\log_2 N \;\text{ การคูณเชิงซ้อน}

N = 2048: จาก 4.2 ล้านเหลือ 11,264 ครั้ง (Cooley–Tukey, N เป็นกำลังของ 2) ผลลัพธ์เหมือน DFT ทุกประการ

ตัวอย่าง D-1 · DFT ด้วยมือ

สุ่มไซน์ $\sin\omega_0 t$ 4 จุดต่อคาบได้ $x = [0, 1, 0, -1]$ หา $X[1]$ แอมพลิจูดและมุม

D1: $e^{-j2\pi n/4} = 1, -j, -1, j$ → $X[1] = 0(1) + 1(-j) + 0(-1) + (-1)(j) = -2j$
D2: $A_1 = 2|{-2j}|/4 = 1$ · มุม −90° (คือ $\cos(\omega t - 90^\circ) = \sin\omega t$ )

แอมพลิจูด 1 มุม −90° ตรงกับไซน์ที่สุ่มมาพอดี แม้มีแค่ 4 จุด เพราะสัญญาณมีความถี่เดียวซึ่งต่ำกว่า $f_s/2$

อเลียสซิง: ฮาร์มอนิกสูงพับลงมาทับฮาร์มอนิกต่ำ กลาง

D5 Harmonic aliasing ฮาร์มอนิกที่เกิน Nyquist

f_s \gt 2f_{max},\qquad h \gt \frac{N}{2} \;\Rightarrow\; \text{ไปปรากฏที่ช่อง } |h - \ell N| \text{ ที่ใกล้ที่สุด}

คลื่นที่มีรอยกระโดด (สี่เหลี่ยม PWM) มีฮาร์มอนิกไม่สิ้นสุด จึงต้องมีตัวกรองกันอเลียสซิงก่อน ADC หรือสุ่มเร็วมาก · เชื่อมกับการสุ่มในสโคป O12 และขบวนอิมพัลส์ C-3

ตัวอย่าง D-2 · สี่เหลี่ยมสุ่มน้อยเกินไป

สุ่มสี่เหลี่ยม ±1 จำนวน 8 จุดและ 64 จุดต่อคาบ DFT ให้ $A_1$ , $A_3$ เท่าไร เทียบค่าจริง $4/\pi = 1.2732$ และ $4/(3\pi) = 0.4244$

8 จุด: $A_1 = 1.3066$ (+2.6 %) · $A_3 = 0.5412$ (+27.5 %) เพราะฮาร์มอนิกที่ 5 พับลงมาที่ 3 ( $8 - 5 = 3$ ), 7 พับลงที่ 1, 9 พับที่ 1, …
64 จุด: $A_1 = 1.2738$ · $A_3 = 0.4260$ ดีขึ้นมาก แต่ $A_{31}$ (ใกล้ Nyquist) ได้ 0.0626 ทั้งที่จริง 0.0411

ฮาร์มอนิกต่ำแม่นเมื่อสุ่มหนาแน่น แต่ช่องที่ใกล้ $f_s/2$ ผิดเสมอ เครื่องวิเคราะห์จึงสุ่ม 256 จุดต่อคาบขึ้นไปและรายงานแค่ถึงลำดับ 50

การรั่วของสเปกตรัม (spectral leakage) และหน้าต่าง กลาง

DFT สมมติว่าข้อมูลในหน้าต่างซ้ำไปเรื่อย ๆ ถ้าหน้าต่างไม่ครบจำนวนเต็มคาบ ต่อปลายกับต้นแล้วเกิดรอยกระโดดที่ไม่มีอยู่จริง (P8) พลังงานจึงกระจายไปทุกช่อง

D6 Leakage and windows การรั่วและหน้าต่าง

\text{สี่เหลี่ยม: } |X[k]| \propto |\mathrm{sinc}(k - k_0)| \qquad \text{Hann: } w[n] = \tfrac12\Big(1 - \cos\frac{2\pi n}{N}\Big),\; \text{ขยายแอมพลิจูด} \times 2

$k_0$ = ตำแหน่งจริงของความถี่ (ไม่เป็นจำนวนเต็ม) · หน้าต่าง Hann ลดการรั่วไปช่องไกล แต่กว้างขึ้น · flat-top อ่านแอมพลิจูดแม่นที่สุดแม้ความถี่อยู่ระหว่างช่อง · ทางที่ดีที่สุดในงานไฟฟ้าคือซิงก์อัตราสุ่มกับความถี่ไฟ (PLL) ให้หน้าต่างครบคาบพอดี

ไซน์ 50.0 Hz (บน) กับ 50.5 Hz (ล่าง) ในหน้าต่าง 200 ms หน้าต่างสี่เหลี่ยม: 50.0 Hz มีช่องเดียว (10) ส่วน 50.5 Hz รั่วไปช่อง 9 และ 11 และไกลออกไป

ตัวอย่าง D-3 · ความถี่ไฟเพี้ยนไป 0.5 Hz

เครื่องวัดสุ่ม 10.24 kHz จำนวน 2048 จุด (200 ms) แต่ไฟจริง 50.5 Hz 230 V หน้าต่างสี่เหลี่ยม อ่านค่าที่ช่อง 9, 10, 11 ได้เท่าไร

หน้าต่างครอบ 10.1 คาบ ไม่ครบจำนวนเต็ม D6
คำนวณ DFT: ช่อง 10 (50 Hz) 225.3 V (−2.0 %) · ช่อง 9 (45 Hz) 19.6 V · ช่อง 11 (55 Hz) 26.0 V ซึ่งไม่มีอยู่จริง
รวมกลุ่มช่อง 9–11 แบบพาร์เซวาล $\sqrt{19.6^2 + 225.3^2 + 26.0^2} = 227.6$ V (−1.0 %)

ความถี่คลาดเพียง 1 % ทำให้ค่าหลักผิด 2 % และเกิดอินเตอร์ฮาร์มอนิกปลอม 8–11 % มาตรฐานจึงบังคับให้ซิงก์กับความถี่ไฟ และใช้การรวมกลุ่ม (grouping) เพื่อลดผลที่เหลือ

มาตรฐานการวัด: IEC 61000-4-7 และ 61000-4-30 กลาง

D7 A typical analyzer setting ตัวอย่างการตั้งค่าเครื่องวิเคราะห์

f_s = 10.24\ \text{kHz},\quad N = 2048 \;\Rightarrow\; T_w = 200\ \text{ms} = 10 \text{ คาบ},\quad \Delta f = 5\ \text{Hz},\quad h \to k = 10h

Nyquist 5.12 kHz สูงกว่าฮาร์มอนิกที่ 50 (2.5 kHz) สองเท่า เผื่อให้ตัวกรองกันอเลียสซิง

D8 IEC 61000-4-7 / 61000-4-30 essentials หลักของมาตรฐานการวัด

T_w = 10 \text{ คาบ (50 Hz) หรือ } 12 \text{ คาบ (60 Hz)} \approx 200\ \text{ms},\qquad G_h = \sqrt{\sum_{k=10h-1}^{10h+1} V_k^2}\;(\text{harmonic subgroup})

อัตราสุ่มซิงก์กับความถี่ไฟ · ค่าราย 200 ms รวมเป็นค่าราย 150 คาบ (≈ 3 s) และ 10 นาที (Class A ของ 61000-4-30) แล้วจึงนำไปเทียบขีดจำกัด เช่น IEEE 519 ใช้เปอร์เซ็นไทล์ของค่าราย 3 s และ 10 นาที (H13) · ช่องระหว่างฮาร์มอนิกคืออินเตอร์ฮาร์มอนิก

ใช้ FFT บนออสซิลโลสโคปวัดฮาร์มอนิกให้ถูก พื้นฐาน

เลือก time/div ให้จอมีหลายคาบ (ความละเอียด $\Delta f = 1/T_{จอ}$ ) แต่ไม่ช้าจนอัตราสุ่มตก (O13)
หน้าต่าง: Hanning สำหรับดูทั่วไป Flat-top เมื่อต้องการแอมพลิจูดแม่น Rectangular เฉพาะเมื่อซิงก์ครบคาบ
หน่วย: ส่วนใหญ่แสดง dBV (RMS) แปลงด้วย D9 แล้วคำนวณ THD เอง หรือใช้ X10
กระแส: ใช้โพรบกระแสหรือ shunt + โพรบดิฟเฟอเรนเชียล ห้ามหนีบกราวด์โพรบที่ไฟบ้าน (ความปลอดภัย)

D9 dBV ↔ V_rms แปลงหน่วยบนจอ FFT

V_{rms} = 10^{L/20}\ \text{V},\qquad L = 20\log_{10}V_{rms}\ \text{dBV}

ตัวอย่าง D-4 · THD จากจอ FFT

FFT ของขาออกเครื่องขยายอ่านได้: หลัก 14 dBV ที่ 3 −2 dBV ที่ 5 −8 dBV หา THD

D9: $V_1 = 10^{0.7} = 5.01$ V · $V_3 = 10^{-0.1} = 0.794$ V · $V_5 = 10^{-0.4} = 0.398$ V
H2: THD $= \sqrt{0.794^2 + 0.398^2}/5.01 = 17.7$ %

THD 17.7 % มีแต่ฮาร์มอนิกคี่ → ขาออกถูกตัดยอดแบบสมมาตร (S7) ลดแอมพลิจูดขาเข้าหรือเพิ่มแรงดันไฟเลี้ยง

คำนวณฮาร์มอนิกบนไมโครคอนโทรลเลอร์: อัลกอริทึม Goertzel ขั้นสูง

เครื่องวัดพลังงานที่ทำเองด้วย ESP32 ไม่จำเป็นต้องใช้ FFT ทั้งก้อน ถ้าต้องการแค่ฮาร์มอนิกที่ 1–15 Goertzel คำนวณทีละช่องด้วยตัวกรอง IIR อันดับสองหนึ่งตัว ทำงานไปพร้อมกับการสุ่มได้ ใช้หน่วยความจำแค่สองตัวแปรต่อช่อง

D10 Goertzel recursion อัลกอริทึม Goertzel

s[n] = x[n] + 2\cos\omega_k\,s[n-1] - s[n-2],\qquad X[k] = s[N-1] - e^{-j\omega_k}s[N-2],\qquad \omega_k = \frac{2\pi k}{N}

เท่ากับ DFT ช่องที่ k ทุกประการ (ต่างกันแค่มุมเฟสคงที่) คูณจริงหนึ่งครั้งต่อตัวอย่าง · หน้าต่างต้องครบคาบ ใช้การวัดความถี่ไฟแล้วปรับอัตราสุ่ม หรือปรับ N ให้ครบ 10 คาบ

// Goertzel: แอมพลิจูด (ค่ายอด) ของช่อง k จาก N ตัวอย่าง — ตรงกับ D2 และ D10
#include <math.h>
float goertzel_amp(const float *x, int N, int k) {
  float w = 2.0f * (float)M_PI * k / N, c = 2.0f * cosf(w);
  float s1 = 0.0f, s2 = 0.0f;
  for (int n = 0; n < N; n++) {
    float s0 = x[n] + c * s1 - s2;
    s2 = s1; s1 = s0;
  }
  float re = s1 - s2 * cosf(w), im = s2 * sinf(w);
  return 2.0f * sqrtf(re * re + im * im) / N;
}

// สุ่ม 128 จุดต่อคาบของ 50 Hz (fs = 6.4 kHz) นาน 10 คาบ → N = 1280, ฮาร์มอนิก h อยู่ช่อง k = 10h
float thd_percent(const float *x) {
  const int N = 1280, M = 10;
  float a1 = goertzel_amp(x, N, M * 1), sum = 0.0f;
  for (int h = 2; h <= 15; h++) { float a = goertzel_amp(x, N, M * h); sum += a * a; }
  return 100.0f * sqrtf(sum) / a1;   // THD จนถึงฮาร์มอนิกที่ 15 (Nyquist = ฮาร์มอนิกที่ 64)
}

ลบค่าเฉลี่ย (DC offset ของ ADC) ออกจาก x ก่อน และกรองกันอเลียสซิงด้วย RC ก่อนเข้า ADC ชิปวัดพลังงานบางรุ่นมีวงจร DFT ในตัว ให้ค่าฮาร์มอนิกโดยตรง

แล็บ: สเปกตรัมแบบเครื่องวัด พื้นฐาน

แล็บในชุดนี้คำนวณสเปกตรัมแบบปิด (อินทิเกรตแต่ละช่วงตรง ๆ) จึงไม่มีอเลียสซิงหรือ leakage ใช้เป็น “คำตอบอ้างอิง” เทียบกับ FFT บนสโคปจริงของคุณ ส่วนแล็บเครื่องกำเนิดสัญญาณมี FFT แบบสโคปจริงที่เกิดอเลียสซิงได้ ลองอเลียสซิง

คำถามที่พบบ่อย

DFT กับอนุกรมฟูเรียร์เกี่ยวกันอย่างไร?

ถ้าสุ่มสัญญาณคาบจำนวน N จุดเท่ากันพอดีในหนึ่งคาบ (หรือจำนวนเต็มคาบ) DFT ของข้อมูลนั้นหารด้วย N คือสัมประสิทธิ์ฟูเรียร์ $c_k$ ตรงตัว ตราบที่สัญญาณไม่มีฮาร์มอนิกเกินครึ่งหนึ่งของความถี่สุ่ม FFT เป็นแค่วิธีคำนวณ DFT ที่เร็ว

ทำไมเครื่องวิเคราะห์กำลังใช้หน้าต่าง 10 คาบ?

IEC 61000-4-7 กำหนดหน้าต่าง 10 คาบที่ 50 Hz (12 คาบที่ 60 Hz) ประมาณ 200 ms ได้ความละเอียด 5 Hz ฮาร์มอนิกที่ h อยู่ช่องที่ 10h และช่องระหว่างนั้นใช้วัดอินเตอร์ฮาร์มอนิก การสุ่มต้องซิงก์กับความถี่ไฟจริงเพื่อไม่ให้เกิด leakage

spectral leakage คืออะไร?

ถ้าหน้าต่างไม่ครอบจำนวนเต็มคาบพอดี DFT มองว่าสัญญาณมีรอยต่อที่ขอบหน้าต่าง พลังงานของฮาร์มอนิกจึงรั่วไปช่องข้างเคียง ค่าแอมพลิจูดผิดและเกิดเส้นปลอม แก้ด้วยการซิงก์อัตราสุ่มกับความถี่ (PLL) ใช้หน้าต่าง Hann หรือ flat-top และการรวมกลุ่มช่อง

Goertzel ดีกว่า FFT เมื่อไร?

เมื่อต้องการแค่ไม่กี่ความถี่ เช่นฮาร์มอนิก 1 ถึง 15 ของไฟ 50 Hz บนไมโครคอนโทรลเลอร์ Goertzel คำนวณแต่ละช่องด้วยลูปเดียว ใช้หน่วยความจำน้อย ทำได้ทีละตัวอย่างระหว่างสุ่ม ส่วน FFT คุ้มกว่าเมื่อต้องการทุกช่อง

อาจารย์มด · ธีรวุฒิ แสวงบุญ
อาจารย์สาขาวิชาวิศวกรรมไฟฟ้า มหาวิทยาลัยราชภัฏธนบุรี สมุทรปราการ · ดูโปรไฟล์

จากอินทิกรัลสู่ผลบวก: DFT คืออนุกรมฟูเรียร์ของข้อมูลที่สุ่มมา กลาง

สุ่มไซน์ sin⁡ω0t\sin\omega_0 tsinω0​t 4 จุดต่อคาบได้ x=[0,1,0,−1]x = [0, 1, 0, -1]x=[0,1,0,−1] หา X[1]X[1]X[1] แอมพลิจูดและมุม

อเลียสซิง: ฮาร์มอนิกสูงพับลงมาทับฮาร์มอนิกต่ำ กลาง

การรั่วของสเปกตรัม (spectral leakage) และหน้าต่าง กลาง

เครื่องวัดสุ่ม 10.24 kHz จำนวน 2048 จุด (200 ms) แต่ไฟจริง 50.5 Hz 230 V หน้าต่างสี่เหลี่ยม อ่านค่าที่ช่อง 9, 10, 11 ได้เท่าไร

มาตรฐานการวัด: IEC 61000-4-7 และ 61000-4-30 กลาง

ใช้ FFT บนออสซิลโลสโคปวัดฮาร์มอนิกให้ถูก พื้นฐาน

FFT ของขาออกเครื่องขยายอ่านได้: หลัก 14 dBV ที่ 3 −2 dBV ที่ 5 −8 dBV หา THD

คำนวณฮาร์มอนิกบนไมโครคอนโทรลเลอร์: อัลกอริทึม Goertzel ขั้นสูง

แล็บ: สเปกตรัมแบบเครื่องวัด พื้นฐาน

คำถามที่พบบ่อย

สุ่มไซน์ $\sin\omega_0 t$ 4 จุดต่อคาบได้ $x = [0, 1, 0, -1]$ หา $X[1]$ แอมพลิจูดและมุม