ตอนนี้เราพร้อมพิจารณาการออกแบบแบบเป็นระบบสำหรับ nonrecursive filters แล้ว วิธีการออกแบบอาศัย Figure 16.1 ซึ่งมีหกส่วน ด้านบนซ้ายเป็นภาพร่างของฟิลเตอร์อุดมคติที่ต้องการ อาจเป็น low-pass, high-pass, band-pass, band-stop, notch filter หรือแม้แต่ differentiator สำหรับฟิลเตอร์ที่ไม่ใช่ differentiator มักต้องการให้ค่าสูงเป็น 0 หรือ 1 ในช่วงต่าง ๆ แต่สำหรับ differentiator ต้องการ เพราะอนุพันธ์ของ eigenfunction เป็น

สรุป: and hence the desired eigenvalues are the coefficient . For a differentiator there is likely to be a cutoff at some frequency because, as you can see, differentiation magnifies, multiplies by ω, and is larger at the high frequencies, which is where the noise usually is, Figure 16.2. See also Figure 15.2.

Figure 16.1—การตัดอนุกรม Fourier แบบอนันต์

สรุป: Figure 16.2—Amplification for differentiation

สัมประสิทธิ์ของอนุกรม Fourier แบบเป็นเชิงรูป (formal Fourier series) ที่สอดคล้องกันคำนวณได้ไม่ยาก เพราะตัวถูกอินทิเกรนด์ในนิพจน์ตรงไปตรงมา (ใช้ integration by parts เมื่อมีอนุพันธ์) สมมติว่าเรานำเสนออนุกรมในรูปของ complex exponentials แล้วสัมประสิทธิ์ของฟิลเตอร์ก็เป็นเพียง Fourier coefficients ของเทอมเอกซ์โพเนนเชียลที่สอดคล้องกัน ด้านขวาบนของ Figure 16.1 แสดงภาพร่างของสัมประสิทธิ์ในเชิงสัญลักษณ์ (แน่นอนว่าเป็นจำนวนเชิงซ้อน)

สรุป: Next, we must truncate the infinite Fourier series to 2N+1 terms (meaning use a rectangular window), shown in Figure 16.1 with the corresponding Fourier representation on the left showing the Gibbs effect.

ขั้นตอนที่สามคือเลือก window เพื่อกำจัด Gibbs effect ที่ร้ายแรงที่สุด สัมประสิทธิ์หลังใส่ window แสดงอยู่ด้านล่างขวา โดยฟิลเตอร์ดิจิทัลสุดท้ายแสดงอยู่ด้านล่างซ้าย ในทางปฏิบัติ ควรกระทำการปัดเศษค่าสัมประสิทธิ์ของฟิลเตอร์ก่อนประเมิน transfer function เพื่อให้เห็นผลของการปัดเศษ

ในวิธีที่สรุปข้างต้น คุณต้องเลือกทั้ง N ซึ่งเป็นจำนวนเทอมที่เก็บไว้ และรูปร่างของ window ถ้าผลที่ได้ไม่เป็นที่พอใจ ก็ต้องเลือกใหม่ นี่เป็นวิธีการออกแบบแบบ "trial and error"

J.F. Kaiser ได้เสนอวิธีการออกแบบที่หาทั้ง N และสมาชิกของตระกูลของ windows ที่ทำงานได้ คุณต้องกำหนดสองอย่างนอกเหนือจากรูปร่าง: ระยะเชิงตั้งฉากที่คุณยอมให้แตกต่างจากไอเดีย (labeled δ) และความกว้างการเปลี่ยนผ่านระหว่าง pass และ stop bands (labeled ΔF) ดู Figure 16.3

Figure 16.3—แบนด์พาสฟิลเตอร์

สำหรับ band-pass filter โดยมี f_p เป็น band-pass และ f_s เป็น band-stop frequencies ลำดับสูตรออกแบบคือ

ถ้า N ใหญ่เกินไป ให้หยุดและพิจารณาการออกแบบใหม่ มิฉะนั้นดำเนินการคำนวณตามลำดับ

สรุป: (this is plotted in Figure 16.4). The original Fourier coefficients for a band-pass filter are given by

สรุป: Figure 16.4—α (A)

สัมประสิทธิ์เหล่านี้จะถูกคูณด้วยน้ำหนักที่สอดคล้องกัน w_k ของ window

สรุป: where

I_0(_x) เป็นฟังก์ชัน Bessel แบบมีส่วนจินตภาพลำดับ 0 ในการคำนวณมันคุณจะต้องการเทอมค่อนข้างน้อย เพราะในตัวส่วนมี (n)!2 ทำให้ลำดับรวมลู่เข้าอย่างรวดเร็ว

สรุป: I_0(_x) is best computed recursively; for a given x, the successive terms of the series are given by

สำหรับ low-pass หรือ high-pass หนึ่งในสองความถี่ f_p หรือ f_s จะมีขอบเขตสูงสุดที่เป็นไปได้ สำหรับ band-stop filter สูตรสัมประสิทธิ์ c_k จะมีการเปลี่ยนแปลงเล็กน้อย

สรุป: Let us examine Kaiser’s window coefficients, the _w_k:

สรุป: As we examine these numbers we see they have, for α > 0, something like the shape of a raised cosine

สรุป: and resemble the von Hann and Hamming windows. There is a “platform” when A > 21. For A > 21 then a = 0, all the w_k are equal to 1, and it is a Lanczos-type window. As _A increases the platform gradually appears. Thus the Kaiser window has properties like many of the more popular ones, and the particular window you use is determined from your specifications via his window rather than by guess or prejudice.

สรุป: How did Kaiser find the formulas? To some extent by trial and error. He first assumed he had a single discontinuity, and he ran a large number of cases on a computer to see both the rise time ΔF and the ripple height δ. With a fair amount of thinking, plus a touch of genius, and noting as a function of A that as A increases we pass from a Lanczos window (A < 21) to a platform of increasing height, 1/I_0 (α). Ideally he wanted a prolate spheroidal function, but he noted they are accurately approximated, for his values, by the _I_0(_x). He plotted the results and approximated the functions. I asked him how he got the exponent 0.4. He replied he tried 0.5 and it was too large, and 0.4, being the next natural choice, seemed to fit very well. It is a good example of using what one knows plus the computer as an experimental tool, even in theoretical research, to get very useful results.

วิธีของ Kaiser จะล้มเหลวบ้างเป็นครั้งคราว เพราะอาจมีมากกว่าหนึ่งขอบ (จริงอยู่ที่มีภาพสะท้อนเชิงสมมาตรของขอบบนส่วนลบของเส้นความถี่) และริพลจากขอบต่าง ๆ อาจบังเอิญรวมกันทำให้ริพลของฟิลเตอร์เกินขนาดที่กำหนด ในกรณีนี้ ซึ่งเกิดขึ้นไม่บ่อยนัก ให้ทำการออกแบบซ้ำด้วยความทนผิดพลาดที่เล็กลง โปรแกรมทั้งหมดสามารถรันได้บนเครื่องคิดเลขโปรแกรมได้แบบพกพาโบราณอย่าง ti–59 น้อยไปกว่าบนพีซีสมัยใหม่

สรุป: We next turn to the finite Fourier series. It is a remarkable fact the Fourier functions are orthogonal, not only over a line segment, but for any discrete set of equally spaced points. Hence the theory will go much the same, except there can be only as many coefficients determined in the Fourier series as there are points. In the case of 2N points, the common case, there is one term of the highest frequency only, the cosine term (the sine term would be identically zero at the sample points). The coefficients are determined as sums of the data points multiplied by the appropriate Fourier functions. The resulting representation will, within roundoff, reproduce the original data.

การคำนวณการขยายจะดูเหมือน 2N เทอมที่แต่ละเทอมมีการคูณและบวก 2N ครั้ง ดังนั้นจึงเป็นการประมาณ (2N)^2 การดำเนินการคูณและบวก แต่โดยทั้ง (1) การบวกและการลบของเทอมที่มีตัวคูณเดียวกันก่อนการคูณ และ (2) การสร้างความถี่ที่สูงขึ้นโดยคูณความถี่ที่ต่ำกว่า จึงทำให้เกิด fast Fourier transform (fft) ที่ต้องการประมาณ N log N การลดความพยายามในการคำนวณนี้ได้เปลี่ยนแปลงพื้นที่ทั้งของวิทยาศาสตร์และวิศวกรรม—สิ่งที่ครั้งหนึ่งเป็นไปไม่ได้ทั้งด้านเวลาและต้นทุน กลับกลายเป็นเรื่องปกติ

สรุป: Now for another story from life. You have all heard about the fast Fourier transform and the Tukey-Cooley paper. It is sometimes called the Tukey-Cooley transform or algorithm. Tukey had suggested to me, sort of, the basic ideas of the fft. I had at that time an ibm Card Programmed Calculator (cpc), and the “butterfly” operation meant it was completely impracticable to do with the equipment I had. Some years later I had an internally programmed ibm 650 and he remarked on it again. All I remembered was it was one of Tukey’s few bad ideas; I completely forgot why it was bad—namely because of the equipment I had at the time. So I did not do the fft, though a book I had already published (1961) shows I knew all the facts necessary, and could have done it easily!

บทสรุป: เมื่อคุณรู้ว่าบางอย่างทำไม่ได้ อย่าลืมเหตุผลหลักว่าทำไม เพื่อว่าเมื่อสถานการณ์เปลี่ยนไปในภายหลัง คุณจะไม่บอกว่า “It can’t be done.” คิดถึงความผิดพลาดของผม! จะโง่ไปกว่านี้อีกได้ไหม? โชคดีที่สำหรับอีโก้ของผม มันเป็นความผิดพลาดที่พบได้บ่อย (และผมเคยทำมากกว่าหนึ่งครั้ง) แต่เนื่องจากความผิดพลาดเกี่ยวกับ fft ทำให้ผมไวต่อเรื่องนี้มากขึ้น ผมยังสังเกตเมื่อคนอื่นทำเช่นกัน—ซึ่งเกิดบ่อยเกินไป! โปรดจำเรื่องความโง่ของผมและสิ่งที่ผมพลาด และอย่าทำความผิดพลาดนั้นเอง เมื่อคุณตัดสินใจว่าสิ่งใดเป็นไปไม่ได้ อย่าพูดทีหลังว่ามันยังเป็นไปไม่ได้โดยไม่ทบทวนทุกรายละเอียดของเหตุผลที่คุณเคยบอกว่ามันทำไม่ได้

ตอนนี้ผมต้องหันไปยังหัวข้ออ่อนไหวเรื่อง power spectra ซึ่งเป็นผลรวมของกำลังสองของสัมประสิทธิ์สองตัวของความถี่หนึ่ง ๆ ในโดเมนจริง หรือกำลังสองของค่าสัมบูรณ์ในสัญกรณ์เชิงซ้อน การพิจารณามันจะทำให้เห็นว่าสิ่งนี้ไม่ขึ้นกับจุดเริ่มต้นของเวลา แต่ขึ้นกับสัญญาณเองเพียงอย่างเดียว ตรงกันข้ามกับการขึ้นกับตำแหน่งของจุดเริ่มต้นของสัมประสิทธิ์ สเปกตรัมมีบทบาทสำคัญในประวัติศาสตร์วิทยาศาสตร์และวิศวกรรม มันเป็น spectral lines ที่เปิดกล่องดำของอะตอมและทำให้ Bohr มองเห็นภายใน กลศาสตร์ควอนตัมใหม่เริ่มรอบปี 1925 ปรับเปลี่ยนเล็กน้อย แต่สเปกตรัมยังคงเป็นกุญแจ เรายังวิเคราะห์กล่องดำเป็นประจำโดยดูสเปกตรัมของอินพุตและเอาต์พุต พร้อมกับความสัมพันธ์ (correlations) เพื่อเข้าใจภายใน—ไม่ได้หมายความว่าจะมีคำตอบเดียวเสมอ แต่โดยทั่วไปเราจะได้เบาะแสเพียงพอที่จะตั้งทฤษฎีใหม่

สรุป: Let us analyze carefully what we do and its implications, because what we do to a great extent controls what we can see. There is usually, in our imaginations at least, a continuous signal. This is usually endless, and we take a sample in time of length 2L. This is the same as multiplying the signal by a Lanczos window, a box car if you prefer. This means the original signal is convolved with the corresponding function of the form (sin x)/x, Figure 16.5—the longer the signal, the narrower the (sin x)/x loops are. Each pure spectral line is smeared out into its (sin x)/x shape.

Figure 16.5—กราฟของ sin x/x

ต่อไปเราทำการสุ่มตัวอย่างที่ช่วงเวลาเท่า ๆ กัน และทุกความถี่ที่สูงกว่าจะถูก alias เข้าเป็นความถี่ที่ต่ำกว่า มันชัดเจนว่าการสลับลำดับของการดำเนินการทั้งสองนี้—การสุ่มแล้วตามด้วยจำกัดช่วง หรือจำกัดช่วงแล้วตามด้วยการสุ่ม—จะให้ผลลัพธ์เหมือนกัน และตามที่ผมกล่าวก่อนหน้านี้ ผมเคยคำนวณพีชคณิตทั้งหมดอย่างละเอียดเพื่อยืนยันว่าทฤษฎีที่ผมคิดต้องเป็นจริงในทางปฏิบัติ

สรุป: Then we use the fft, which is only a cute, accurate way of getting the coefficients of a finite Fourier series. But when we assume the finite Fourier series representation _we are making the function periodic_—and the period is exactly the sampling interval size times the number of samples we take! This period has generally nothing to do with the periods in the original signal. _We force all non-harmonic frequencies into harmonic ones_—we force a continuous spectrum to be a line spectrum! This forcing is not a local effect, but as you can easily compute, a non-harmonic frequency goes into all the other frequencies, most strongly into the adjacent ones of course, but nontrivially into more remote frequencies.

ผมได้ข้ามเทคนิคสถิติมาตรฐานในการตัดค่าเฉลี่ยออก ไม่ว่าจะเพื่อความสะดวกหรือเหตุผลการคาลิเบรต ซึ่งจะลดพลังงานที่ความถี่ 0 ในสเปกตรัมให้เป็น 0 และสร้างความไม่ต่อเนื่องอย่างสำคัญในสเปกตรัม หากคุณใช้ window ทีหลัง มันเพียงแค่เบลอมันไปรอบ ๆ ความถี่ข้างเคียง ในการประมวลผลข้อมูลสำหรับ Tukey ผมมักตัดเส้นแนวโน้มเชิงเส้นและแม้แต่พาราโบลาแนวโน้มจากข้อมูลเที่ยวบินของเครื่องบินหรือขีปนาวุธ แล้ววิเคราะห์เศษที่เหลือ แต่สเปกตรัมของผลบวกของสองสัญญาณไม่ใช่ผลบวกของสเปกตรัม—ไม่ใช่เลย! เมื่อคุณบวกสองฟังก์ชัน ความถี่ต่าง ๆ จะถูกบวก เชิงพีชคณิต และพวกมันอาจเสริมกันหรือลดทอนกัน และดังนั้นให้ผลลัพธ์ที่ผิดพลาดโดยสิ้นเชิง! ไม่มีใครที่ผมรู้จักมีคำตอบที่สมเหตุสมผลต่อข้อโต้แย้งของผมที่นี่—เรายังทำมันบางส่วนเพราะเราไม่รู้จะทำอย่างอื่น—แต่เส้นแนวโน้มมีความไม่ต่อเนื่องขนาดใหญ่ที่ปลาย (จำไว้ว่าพวกเราสมมติว่าฟังก์ชันทั้งหมดเป็นแบบ periodic) และดังนั้นสัมประสิทธิ์ของมันจึงลดลงเป็น 1/k ซึ่งไม่ลดอย่างรวดเร็วเลย

สรุป: Let us turn to theory. Every spectrum of real noise falls off reasonably rapidly as you go to infinite frequencies, or else it would have infinite energy, Figure 16.6. But the sampling process aliases the higher frequencies in lower ones, and the folding as indicated tends to produce a flat spectrum—remember, the frequencies when aliased are algebraically added. Hence we tend to see a flat spectrum for noise, and if it is flat then we call it white noise. The signal, usually, is mainly in the lower frequencies. This is true for several reasons, including “over-sampling” (sampling more often than is required from the Nyquist theorem), and means we can get some averaging to reduce the instrumental errors. Thus the typical spectrum will look as shown in Figure 16.6. Hence the prevalence of low-pass filters to remove the noise. No linear method can separate the signal from the noise at the same frequencies, but those beyond the signal can be so removed by a low-pass filter. Therefore, when we “over-sample” we have a chance to remove more of the noise by a low-pass filter.

Figure 16.6—สเปกตรัมทั่วไป

จำไว้ว่ามีความเข้าใจโดยนัยว่าเรากำลังประมวลผลระบบเชิงเส้น การวิเคราะห์ Fourier ตลาดหุ้นเก่าที่เผยว่ามีเพียง white noise ถูกตีความว่าไม่มีวิธีทำนายราคาหุ้นในอนาคต—และนี่ถูกต้องเฉพาะเมื่อคุณตั้งใจใช้ predictor เชิงเส้นอย่างง่าย มันไม่ได้บอกอะไรเกี่ยวกับการใช้ predictor ที่ไม่เชิงเส้น อย่างไรก็ตาม อีกครั้งเป็นการตีความผลที่กว้างขวางผิดพลาดเพราะขาดความเข้าใจพื้นฐานของเครื่องมือคณิตศาสตร์ และรู้เพียงแค่เครื่องมือเอง ความรู้เล็กน้อยเป็นสิ่งอันตราย—โดยเฉพาะถ้าคุณขาดพื้นฐาน!

ผมบอกอย่างระมัดระวังในการบรรยายเปิดเกี่ยวกับ digital filters ว่า ณ เวลานั้นผม คิดว่า ผมไม่รู้อะไรเกี่ยวกับพวกมัน สิ่งที่ผมไม่รู้คือเพราะผมไม่รู้เรื่องการออกแบบ recursive digital filter ตอนนั้นผมได้สร้างมันขึ้นมาโดยปริยายเมื่อผมตรวจทฤษฎีของ predictor-corrector methods ในการแก้สมการเชิงอนุพันธ์ธรรมดาโดยตัวเลขอย่างละเอียด corrector นั้นเป็นฟิลเตอร์ดิจิทัลเชิง recursive จริง ๆ

ในขณะที่ศึกษาวิธีอินทิเกรตระบบสมการเชิงอนุพันธ์ธรรมดาด้วยวิธีตัวเลข ผมไม่ได้ถูกรบกวนโดยความคิดเดิม ๆ เกี่ยวกับ digital filters และผมตระหนักไว้ว่าการป้อนข้อมูลที่มีขอบเขต (bounded input) ตามคำพูดของผู้เชี่ยวชาญด้านฟิลเตอร์ อาจทำให้เกิดผลลัพธ์ที่ไม่มีขอบเขต (unbounded output) ถ้าคุณกำลังอินทิเกรต—ซึ่งพวกเขาเรียกว่า unstable —แต่ชัดเจนว่ามันเป็นสิ่งที่คุณต้องการเมื่อต้องอินทิเกรต; แม้แต่ค่าคงที่ก็จะทำให้เกิดการเติบโตเชิงเส้นในผลลัพธ์ ต่อมาเมื่อผมต้องอินทิเกรตเส้นทางลงสู่พื้นผิวของดวงจันทร์ ที่ซึ่งไม่มีอากาศ จึงไม่มี drag ไม่มีเทอมอนุพันธ์อันดับหนึ่งชัดเจนในสมการ และผมต้องการนำข้อได้เปรียบนี้มาใช้ด้วยสูตรการอินทิเกรตทางตัวเลขที่เหมาะสม ผมพบว่าต้องยอมให้มีการเติบโตของข้อผิดพลาดเป็นกำลังสอง; ข้อผิดพลาดการปัดเศษเล็ก ๆ ในการคำนวณความเร่งจะไม่ถูกแก้ไขและจะนำไปสู่ข้อผิดพลาดเชิงกำลังสองในตำแหน่ง—ข้อผิดพลาดในความเร่งทำให้ตำแหน่งเติบโตแบบกำลังสอง นั่นคือธรรมชาติของปัญหา แตกต่างจากบนโลก ที่ drag ของอากาศให้ feedback แก้ค่าความเร่งที่ผิดพลาด และจึงแก้ไขข้อผิดพลาดในตำแหน่งบางส่วน

ดังนั้นจนถึงปัจจุบันผมมีทัศนคติว่า stability in digital filters หมายถึง "ไม่เติบโตแบบเอ็กซ์โพเนนเชียล" จาก bounded inputs แต่ยอมให้การเติบโตเป็นพหุนามได้ และนี่ไม่ใช่มาตรฐานเกณฑ์ความมั่นคงที่มาจากฟิลเตอร์แอนะล็อกแบบคลาสสิก ซึ่งถ้าหากไม่ถูกจำกัดคุณอาจทำให้สิ่งต่าง ๆ ร้อนละลายได้—และอย่างไรก็ตาม พวกเขาไม่ได้คิดอย่างจริงจังเกี่ยวกับการอินทิเกรตว่าเป็นกระบวนการฟิลเตอร์

เราจะหยิบหัวข้อสำคัญนี้ของ recursive filters ซึ่งจำเป็นสำหรับการอินทิเกรต ในบทถัดไป