W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】六

2016-05-10 懸鉤子

人耳聽旋律，縱然不能分辨音高時值，卻知音符總有個先來後到。傳說莫札特絕對音感之好，即使是初聞的曲調，也能立刻譜出！！為求一次得到『兩個領域』── 時間、頻率 ── 的好處，於是人們設想了

短時距傅立葉變換

短時距傅立葉變換是傅立葉變換的一種變形，為時頻分析中其中一個重要的工具。

與傅立葉轉換在概念上的區別

將訊號做傅立葉變換後得到的結果，並不能給予關於訊號頻率隨時間改變的任何資訊。以下的例子作為說明：

x(t)=\begin{cases} \cos(440 \pi t); & t < 0.5 \\ \cos(660 \pi t); & 0.5 \leq t < 1 \\ \cos(524 \pi t); & t \ge 1 \end{cases}

傅立葉變換後的頻譜和短時距傅立葉轉換後的結果如下：

傅立葉轉換後, 橫軸為頻率(赫茲)

短時距傅立葉轉換, 橫軸為時間(秒)，縱軸為頻率(赫茲)

由上圖可發現，傅立葉轉換只提供了有哪些頻率成份的資訊，卻沒有提供時間資訊；而短時傅立葉轉換則清楚的提供這兩種資訊。這種時頻分析的方法有利於頻率會隨著時間改變的訊號，如音樂訊號和語音訊號等分析。

定義

數學定義

簡單來說，在連續時間的例子，一個函數可以先乘上僅在一段時間不為零的窗函數再進行一維的傅立葉變換。再將這個窗函數沿著時間軸挪移，所得到一系列的傅立葉變換結果排開則成為二維表象。數學上，這樣的操作可寫為：

X(t, f) = \int_{-\infty}^{\infty} w(t-\tau)x(\tau) e^{-j 2 \pi f \tau} \, d\tau

另外也可用角頻率來表示：

X(t, \omega) = \int_{-\infty}^{\infty} w(t-\tau)x(\tau) e^{-j \omega \tau} \, d\tau

其中 $w(t)$ 是窗函數，窗函數種類有很多種，會在稍後再做仔細討論。 $x(t)$ 是待變換的訊號。 $X(t,\omega)$ 是 $w(t-\tau)x(\tau)$ 的傅立葉變換。隨著 $t$ 的改變，窗函數在時間軸上會有位移。經 $w(t-\tau)x(\tau)$ 後，訊號只留下了窗函數截取的部分做最後的傅立葉轉換。

而反短時距傅立葉轉換，其數學類似傅立葉轉換，但須消除窗函數的作用：

x(t)=w(t_1-t)^{-1} \int_{-\infty}^{\infty} X(t_1, f) e^{j 2 \pi f t}\, df ; w(t_1-t)\ne 0

窗函數

窗函數通常滿足下列特性:

$w(t) = w(-t) \,$ ，即為偶函數。
$max(w(t))=w(0) \,$ ，即窗函數的中央通常是最大值的位置。
$w(t_1)\ge w(t_2), |t_2| \ge |t_1|$ ，即窗函數的值由中央開始向兩側單調遞減。
$w(t)\cong 0 , |t|\to \infty$ ，即窗函數的值向兩側遞減為零。

常見的窗函數有：方形、三角形、高斯函數等，而短時距傅立葉轉換也因窗函數的不同而有不同的名稱。而加伯轉換，即為窗函數是高斯函數的短時距傅立葉轉換，通常沒有特別說明的短時距傅立葉轉換，即為加伯轉換。

……

頻譜(Spectrogram)

Spectrogram即短時傅立葉轉換後結果的絕對值平方，兩者本質上是相同的，在文獻上也常出現spectrogram這個名詞。

SP_x(t,f) = |X(t,f)|^2 = | \int_{-\infty}^{\infty} w(t-\tau)x(\tau) e^{-j 2 \pi f \tau} \, d\tau |^2

───

誰知卻與『測不準原理』不期而遇︰

相對論將『觀察者』帶入物理，改變了『量測』的基本觀念。雖然無限精準的『測量』即使作不到，尚且還可以想像。但是量子力學把『量測』的『測不準』原理放進物理，就是說連想像『粒子』的『軌跡』在原理上都不『允許』！！量子力學是使用著『運算子』operator 的語言來描寫微小粒子之『事件概率』的『波函數』。那『測不準原理』是什麼呢？所謂『經典物理學』classic physics 對一個『物體』運動軌跡的描述是由它的『位置』和『速度』或說『動量』所確定的，一九二七年德國的維爾納‧海森堡 Werner Heisenberg 卻講任何『量子系統』之『量測』必為如下的關係式所制約︰

$\Delta x \Delta p \ge \hbar$

$\Delta t \Delta E \ge \hbar$

這並不是因為觀察者的量測，影響了系統── 比方說用粗大的溫度計量一小杯水的溫度 ──所導致的『觀察者效應』，而是宇宙的本質如此。所以即使是想像一個箱子裡的『電子軌跡』都沒有『旨趣』，你不量測它想說它是『波』或者是『粒子』之象純屬『無謂』。這引發一些物理學大方家不滿，認為量子力學根本尚未『完備』。就像發展完成量子力學『波動方程式』的埃爾溫‧薛丁格，他卻也是提出一個稱之為『薛丁格貓 』之想像實驗的人，用以表達目前量子力學之『哥本哈根詮釋』所必須思考的嚴峻性矛盾問題︰

薛丁格是如此描述這個實驗的：

實驗者甚至可以設置出相當荒謬的案例來。把一隻貓關在一個封閉的鐵容器裡面，並且裝置以下儀器（注意必須保固這儀器不被容器中的貓直接干擾）：在一台蓋革計數器內置入極少量放射性物質，由於物質的數量極少，在一小時內，這個放射性物質至少有一個原子衰變的機率為 50%，它沒有任何原子衰變的機率也同樣為50%；假若衰變事件發生了，則蓋革計數管會放電，通過繼電器啟動一個榔頭，榔頭會打破裝有氰化氫的燒瓶。經過一小時以後，假若沒有發生衰變事件，則貓仍舊存活；否則發生衰變，這套機構被觸發，氰化氫揮發，導致貓隨即死亡。用以描述整個事件的波函數竟然表達出了活貓與死貓各半糾合在一起的狀態。
類似這典型案例的眾多案例裏，原本只局限於原子領域的不明確性被以一種巧妙的機制變為宏觀不明確性，只有通過打開這個箱子來直接觀察才能解除這樣的不明確性。它使得我們難以如此天真地接受採用這種籠統的模型來正確代表實體的量子特性。就其本身的意義而言，它不會蘊含任何不清楚或矛盾的涵義。但是，在一張搖晃或失焦的圖片與雲堆霧層的快照之間，實則有很大的不同之處。

不僅如此，在量子系統中，假使兩個粒子在經過短時間彼此間耦合之後，儘管將這兩個粒子分隔很遠的一段距離，量測其中任何一個粒子，會不可避免地影響到另外一個粒子的度量性質，彷彿有隔空的傳心術一般，這種關聯現象稱之為『量子糾纏 』quantum entanglement 。

─── 摘自《測不準原理》

那個 $\Delta t \Delta E \ge \hbar$ 限制了 STFT 之『解析度』︰

Short-time Fourier transform

The short-time Fourier transform (STFT), or alternatively short-term Fourier transform, is a Fourier-related transform used to determine the sinusoidal frequency and phase content of local sections of a signal as it changes over time.^[1] In practice, the procedure for computing STFTs is to divide a longer time signal into shorter segments of equal length and then compute the Fourier transform separately on each shorter segment. This reveals the Fourier spectrum on each shorter segment. One then usually plots the changing spectra as a function of time.

Impact_STFT

Example of short time Fourier transforms used to determine time of impact from audio signal.

……

Resolution issues

Further information: Gabor limit

One of the pitfalls of the STFT is that it has a fixed resolution. The width of the windowing function relates to how the signal is represented—it determines whether there is good frequency resolution (frequency components close together can be separated) or good time resolution (the time at which frequencies change). A wide window gives better frequency resolution but poor time resolution. A narrower window gives good time resolution but poor frequency resolution. These are called narrowband and wideband transforms, respectively.

Comparison of STFT resolution. Left has better time resolution, and right has better frequency resolution.

This is one of the reasons for the creation of the wavelet transform and multiresolution analysis, which can give good time resolution for high-frequency events and good frequency resolution for low-frequency events, the combination best suited for many real signals.

This property is related to the Heisenberg uncertainty principle, but not directly – see Gabor limit for discussion. The product of the standard deviation in time and frequency is limited. The boundary of the uncertainty principle (best simultaneous resolution of both) is reached with a Gaussian window function, as the Gaussian minimizes the Fourier uncertainty principle. This is called the Gabor transform (and with modifications for multiresolution becomes the Morlet wavelet transform).

One can consider the STFT for varying window size as a two-dimensional domain (time and frequency), as illustrated in the example below, which can be calculated by varying the window size. However, this is no longer a strictly time–frequency representation – the kernel is not constant over the entire signal.

───

Fourier transform on Euclidean space

The Fourier transform can be defined in any arbitrary number of dimensions n. As with the one-dimensional case, there are many conventions. For an integrable function f(x), this article takes the definition:

\hat{f}(\boldsymbol{\xi}) = \mathcal{F}(f)(\boldsymbol{\xi}) = \int_{\R^n} f(\mathbf{x}) e^{-2\pi i \mathbf{x}\cdot\boldsymbol{\xi}} \, d\mathbf{x}

where x and ξ are n-dimensional vectors, and x · ξ is the dot product of the vectors. The dot product is sometimes written as $\left\langle \mathbf x, \boldsymbol \xi \right\rangle$ .

All of the basic properties listed above hold for the n-dimensional Fourier transform, as do Plancherel’s and Parseval’s theorem. When the function is integrable, the Fourier transform is still uniformly continuous and the Riemann–Lebesgue lemma holds. (Stein & Weiss 1971)

Uncertainty principle

For more details on this topic, see Gabor limit.

Generally speaking, the more concentrated f(x) is, the more spread out its Fourier transform $\hat f(\xi)$ must be. In particular, the scaling property of the Fourier transform may be seen as saying: if we “squeeze” a function in x, its Fourier transform “stretches out” in ξ. It is not possible to arbitrarily concentrate both a function and its Fourier transform.

The trade-off between the compaction of a function and its Fourier transform can be formalized in the form of an uncertainty principle by viewing a function and its Fourier transform as conjugate variables with respect to the symplectic form on the time–frequency domain: from the point of view of the linear canonical transformation, the Fourier transform is rotation by 90° in the time–frequency domain, and preserves the symplectic form.

Suppose f(x) is an integrable and square-integrable function. Without loss of generality, assume that f(x) is normalized:

\int_{-\infty}^\infty |f(x)|^2 \,dx=1.

It follows from the Plancherel theorem that $\hat f(\xi)$ is also normalized.

The spread around x = 0 may be measured by the dispersion about zero (Pinsky 2002, p. 131) defined by

D_0(f)=\int_{-\infty}^\infty x^2|f(x)|^2\,dx.

In probability terms, this is the second moment of |f(x)|² about zero.

The Uncertainty principle states that, if f(x) is absolutely continuous and the functions x·f(x) and f′(x) are square integrable, then

D_0(f)D_0(\hat{f}) \geq \frac{1}{16\pi^2}

(Pinsky 2002).

The equality is attained only in the case $f(x)=C_1 \, e^{-\pi x^2/\sigma^2}$ (hence $\hat{f}(\xi)= \sigma C_1 \, e^{-\pi\sigma^2\xi^2}$ ) where σ > 0 is arbitrary and $C_1 = \sqrt[4]{2} / \sqrt{\sigma}$ so that f is L²–normalized (Pinsky 2002). In other words, where f is a (normalized) Gaussian function with variance σ², centered at zero, and its Fourier transform is a Gaussian function with variance σ⁻².

In fact, this inequality implies that:

\left(\int_{-\infty}^\infty (x-x_0)^2|f(x)|^2\,dx\right)\left(\int_{-\infty}^\infty(\xi-\xi_0)^2|\hat{f}(\xi)|^2\,d\xi\right)\geq \frac{1}{16\pi^2}

for any x₀, ξ₀ ∈ R (Stein & Shakarchi 2003, p. 158).

In quantum mechanics, the momentum and position wave functions are Fourier transform pairs, to within a factor of Planck’s constant. With this constant properly taken into account, the inequality above becomes the statement of the Heisenberg uncertainty principle (Stein & Shakarchi 2003, p. 158).

A stronger uncertainty principle is the Hirschman uncertainty principle, which is expressed as:

H(|f|^2)+H(|\hat{f}|^2)\ge \log(e/2)

where H(p) is the differential entropy of the probability density function p(x):

H(p) = -\int_{-\infty}^\infty p(x)\log(p(x)) \, dx

where the logarithms may be in any base that is consistent. The equality is attained for a Gaussian, as in the previous case.

───

由於 numpy 中並沒有直接計算 STFT 的函式，有興趣的讀者或可以安裝 scipy

scipy.signal.spectrogram

scipy.signal.spectrogram(x, fs=1.0, window=(‘tukey’, 0.25), nperseg=256, noverlap=None, nfft=None, detrend=’constant’, return_onesided=True, scaling=’density’, axis=-1, mode=’psd’)¶

Compute a spectrogram with consecutive Fourier transforms.

Spectrograms can be used as a way of visualizing the change of a nonstationary signal’s frequency content over time.

Parameters:

Parameters:	x : array_like Time series of measurement values fs : float, optional Sampling frequency of the x time series. Defaults to 1.0. window : str or tuple or array_like, optional Desired window to use. See `get_window` for a list of windows and required parameters. If window is array_like it will be used directly as the window and its length will be used for nperseg. Defaults to a Tukey window with shape parameter of 0.25. nperseg : int, optional Length of each segment. Defaults to 256. noverlap : int, optional Number of points to overlap between segments. If None, `noverlap = nperseg // 8`. Defaults to None. nfft : int, optional Length of the FFT used, if a zero padded FFT is desired. If None, the FFT length is nperseg. Defaults to None. detrend : str or function or False, optional Specifies how to detrend each segment. If `detrend` is a string, it is passed as the `type` argument to `detrend`. If it is a function, it takes a segment and returns a detrended segment. If `detrend` is False, no detrending is done. Defaults to ‘constant’. return_onesided : bool, optional If True, return a one-sided spectrum for real data. If False return a two-sided spectrum. Note that for complex data, a two-sided spectrum is always returned. scaling : { ‘density’, ‘spectrum’ }, optional Selects between computing the power spectral density (‘density’) where Pxx has units of V*2/Hz and computing the power spectrum (‘spectrum’) where Pxx* has units of V2, if x is measured in V and fs is measured in Hz. Defaults to ‘density’ axis : int, optional Axis along which the spectrogram is computed; the default is over the last axis (i.e. `axis=-1`). mode** : str, optional Defines what kind of return values are expected. Options are [‘psd’, ‘complex’, ‘magnitude’, ‘angle’, ‘phase’].
Returns:	f : ndarray Array of sample frequencies. t : ndarray Array of segment times. Sxx : ndarray Spectrogram of x. By default, the last axis of Sxx corresponds to the segment times.

x : array_like

Time series of measurement values

fs : float, optional

Sampling frequency of the x time series. Defaults to 1.0.

window : str or tuple or array_like, optional

Desired window to use. See get_window for a list of windows and required parameters. If window is array_like it will be used directly as the window and its length will be used for nperseg. Defaults to a Tukey window with shape parameter of 0.25.

nperseg : int, optional

Length of each segment. Defaults to 256.

noverlap : int, optional

Number of points to overlap between segments. If None, noverlap = nperseg // 8. Defaults to None.

nfft : int, optional

Length of the FFT used, if a zero padded FFT is desired. If None, the FFT length is nperseg. Defaults to None.

detrend : str or function or False, optional

Specifies how to detrend each segment. If detrend is a string, it is passed as the type argument to detrend. If it is a function, it takes a segment and returns a detrended segment. If detrend is False, no detrending is done. Defaults to ‘constant’.

return_onesided : bool, optional

If True, return a one-sided spectrum for real data. If False return a two-sided spectrum. Note that for complex data, a two-sided spectrum is always returned.

scaling : { ‘density’, ‘spectrum’ }, optional

Selects between computing the power spectral density (‘density’) where Pxx has units of V**2/Hz and computing the power spectrum (‘spectrum’) where Pxx has units of V**2, if x is measured in V and fs is measured in Hz. Defaults to ‘density’

axis : int, optional

Axis along which the spectrogram is computed; the default is over the last axis (i.e. axis=-1).

mode : str, optional

Defines what kind of return values are expected. Options are [‘psd’, ‘complex’, ‘magnitude’, ‘angle’, ‘phase’].

Returns:

f : ndarray

Array of sample frequencies.

t : ndarray

Array of segment times.

Sxx : ndarray

Spectrogram of x. By default, the last axis of Sxx corresponds to the segment times.

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】五

2016-05-09 懸鉤子

竹窗‧白居易

常愛輞川寺，竹窗東北廊。
一别十馀載，見竹未曾忘。
今春二月初，蔔居在新昌。
未暇作廄庫，且先營一堂。
開窗不糊紙，種竹不依行。
意取北檐下，窗與竹相當。
繞屋聲淅淅，逼人色蒼蒼。
煙通杳靄氣，月透玲瓏光。
是時三伏天，天氣熱如湯。
獨此竹窗下，朝回解衣裳。
輕紗一幅巾，小簟六尺床。
無客盡日靜，有風終夜涼。
乃知前古人，言事頗諳詳。
清風北窗臥，可以傲羲皇。

柏拉圖說︰眼睛是靈魂之窗。因此靈魂無法一眼盡收宇宙萬有，很容易墜入『洞見』的耶？？

【柏拉圖《理想國》第七卷『洞穴寓言』節錄】

蘇：接下來讓我們把受過教育的人與沒受過教育的人的本質比作下述情形。讓我們想像一個洞穴式的地下室，它有一長長通道通向外面，可讓和洞穴一樣寬的一路亮光照進來。有一些人從小就住在這洞穴裡，頭頸和腿腳都綁著，不能走動也不能轉頭，只能向前看著洞穴後壁。讓我們再想像在他們背後遠處高些的地方有東西燃燒著發出火光。在火光和這些被囚禁者之間，在洞外上面有一條路。沿著路邊已築有一帶矮牆。矮牆的作用象傀儡戲演員在自己和觀眾之間設的一道屏障，他們把木偶舉到屏障上頭去表演。
格：我看見了。
蘇：接下來讓我們想像有一些人拿著各種器物舉過牆頭，從牆後面走過，有的還舉著用木料、石料或其它材料製作的假人和假獸。而這些過路人，你可以料到有的在說話，有的不在說話。
格：你說的是一個奇特的比喻和一些奇特的囚徒。
蘇：不，他們是一些和我們一樣的人。你且說說看，你認為這些囚徒除了火光投射到他們對面洞壁上的陰影而外，他們還能看到自己的或同伴們的什麼呢？
格：如果他們一輩子頭頸被限制了不能轉動，他們又怎樣能看到別的什麼呢？
蘇：那麼，後面路上人舉著過去的東西，除了它們的陰影而外，囚徒們能看到它們別的什麼嗎？
格：當然不能。
蘇：那麼，如果囚徒們能彼此交談，你不認為，他們會斷定，他們在講自己所看到的陰影時是在講真物本身嗎？
格：必定如此。
蘇：又，如果一個過路人發出聲音，引起囚徒對面洞壁的回聲，你不認為，囚徒們會斷定，這是他們對面洞壁上移動的陰影發出的嗎？
格：他們一定會這樣斷定的。
蘇：因此無疑，這種人不會想到，上述事物除陰影而外還有什麼別的實在。
格：無疑的。
蘇：那麼，請設想一下，如果他們被解除禁錮，矯正迷誤，你認為這時他們會怎樣呢？如果真的發生如下的事情：其中有一人被解除了桎梏，被迫突然站了起來，轉頭環視，走動，抬頭看望火光，你認為這時他會怎樣呢？他在做這些動作時會感覺痛苦的，並且，由於眼花潦亂，他無法看見那些他原來只看見其陰影的實物。如果有人告訴他，說他過去慣常看到的全然是虛假，如今他由於被扭向了比較真實的器物，比較地接近了實在，所見比較真實了，你認為他聽了這話會說些什麼呢？如果再有人把牆頭上過去的每一器物指給他看，並且逼他說出那是些什麼，你不認為，這時他會不知說什麼是好，並且認為他過去所看到的陰影比現在所看到的實物更真實嗎？
格：更真實得多呀！
蘇：如果他被迫看火光本身，他的眼睛會感到痛苦，他會轉身走開，仍舊逃向那些他能夠看清而且確實認為比人家所指示的實物還更清楚更實在的影像的。不是嗎？
格：會這樣的。
蘇：再說，如果有人硬拉他走上一條陡峭崎嶇的坡道，直到把他拉出洞穴見到了外面的陽光，不讓他中途退回去，他會覺得這樣被強迫著走很痛苦，並且感到惱火；當他來到陽光下時，他會覺得眼前金星亂蹦金蛇亂串，以致無法看見任何一個現在被稱為真實的事物的。你不認為會這樣嗎？
格：噢，的確不是一下子就能看得見的。
蘇：因此我認為，要他能在洞穴外面的高處看得見東西，大概需要有一個逐漸習慣的過程。首先大概看陰影是最容易，其次要數看人和其他東西在水中的倒影容易，再次是看東西本身；經過這些之後他大概會覺得在夜裡觀察天象和天空本身，看月光和星光，比白天看太陽和太陽光容易。
格：當然囉。
蘇：這樣一來，我認為，他大概終於就能直接觀看太陽本身，看見他的真相了，就可以不必通過水中的倒影或影像，或任何其他媒介中顯示出的影像看它了，就可以在它本來的地方就其本身看見其本相了。
格：這是一定的。
蘇：接著他大概對此已經可以得出結論了：造成四季交替和年歲週期，主宰可見世界一切事物的正是這個太陽，它也就是他們過去通過某種曲折看見的所有那些事物的原因。
格：顯然，他大概會接著得出這樣的結論。
蘇：如果他回想自己當初的穴居、那個時候的智力水平，以及禁錮中的夥伴們，你不認為，他會慶幸自己的這一變遷，而替夥伴們遺憾嗎？
格：確實會的。
蘇：如果囚徒們之間曾有過某種選舉，也有人在其中贏得過尊榮，而那些敏於辨別而且最能記住過往影像的慣常次序，因而最能預言後面還有什麼影像會跟上來的人還得到過獎勵，你認為這個既已解放了的人他會再熱衷於這種獎賞嗎？對那些受到囚徒們尊重並成了他們領袖的人，他會心懷嫉妒，和他們爭奪那裡的權力地位嗎？或者，還是會像荷馬所說的那樣，他寧願活在人世上做一個窮人的奴隸，受苦受難，也不願和囚徒們有共同意見，再過他們那種生活呢？
格：我想，他會寧願忍受任何苦楚也不願再過囚徒生活的。
蘇：如果他又回到地穴中坐在他原來的位置上，你認為會怎麼樣呢？他由於突然地離開陽光走進地穴，他的眼睛不會因黑暗而變得什麼也看不見嗎？
格：一定是這樣的。
蘇：這時他的視力還很模糊，還沒來得及習慣於黑暗——再習慣於黑暗所需的時間也不會是很短的。如果有人趁這時就要他和那些始終禁錮在地穴中的人們較量一下「評價影像」，他不會遭到笑話嗎？人家不會說他到上面去走了一趟，回來眼睛就壞了，不會說甚至連起一個往上去的念頭都是不值得的嗎？要是把那個打算釋放他們並把他們帶到上面去的人逮住殺掉是可以的話，他們不會殺掉他嗎？
格：他們一定會的。
蘇：親愛的格勞孔，現在我們必須把這個比喻整個兒地應用到前面講過的事情上去，把地穴囚室比喻可見世界，把火光比喻太陽的能力。如果你把從地穴到上面世界並在上面看見東西的上升過程和靈魂上升到可知世界的上升過程聯想起來，你就領會對了我的這一解釋了，既然你急於要聽我的解釋。至於這一解釋本身是不是對，這是只有神知道的。但是無論如何，我覺得，在可知世界中最後看見的，而且是要花很大的努力才能最後看見的東西乃是善的理念。我們一旦看見了它，就必定能得出下述結論：它的確就是一切事物中一切正確者和美者的原因，就是可見世界中創造光和光源者，在可理知世界中它本身就是真理和理性的決定性源泉；任何人凡能在私人生活或公共生活中行事合乎理性的，必定是看見了善的理念的。
格：就我所能瞭解的而言，我都同意。
蘇：那麼來吧，你也來同意我下述的看法吧，而且在看到下述情形時別感到奇怪吧：那些已達到這一高度的人不願意做那些瑣碎俗事，他們的心靈永遠渴望逗留在高處的真實之境。如果我們的比喻是合適的話，這種情形應該是不奇怪的。
格：是不足為怪的。
蘇：再說，如果有人從神聖的觀察再回到人事；他在還看不見東西還沒有變得足夠地習慣於黑暗環境時，就被迫在法庭上或其它什麼地方同人家爭訟關於正義的影子或產生影子的偶像，辯論從未見過正義本身的人頭腦裡關於正義的觀念。如果他在這樣做時
顯得樣子很難看舉止極可笑，你認為值得奇怪嗎？
格：一點也不值得奇怪。
蘇：但是，凡有頭腦的人都會記得，眼睛有性質不同的兩種迷盲，它們是由兩種相應的原因引起的：一是由亮處到了暗處，另一是由暗處到了亮處。凡有頭腦的人也都會相信，靈魂也能出現同樣的情況。他在看到某個靈魂發生迷盲不能看清事物時，不會不加思索就予以嘲笑的，他會考察一下，靈魂的視覺是因為離開了較光明的生活被不習慣的黑暗迷誤了的呢，還是由於離開了無知的黑暗進入了比較光明的世界，較大的亮光使它失去了視覺的呢？於是他會認為一種經驗與生活道路是幸福的，另一種經驗與生活道路是可憐的；如果他想笑一笑的話，那麼從下面到上面去的那一種是不及從上面的亮處到下面來的這一種可笑的。
格：你說的非常有道理。
蘇：如果這是正確的，那麼關於這些事，我們就必須有如下的看法：教育實際上並不像某些人在自己的職業中所宣稱的那樣。他們宣稱，他們能把靈魂裡原來沒有的知識灌輸到靈魂裡去，好像他們能把視力放進瞎子的眼睛裡去似的。
格：他們確曾有過這種說法。
蘇：…

………

─── 摘自《知識是什麼？》

詩人剪裁天地，衷情自然之美，所以講︰

一别十馀載，見竹未曾忘。

當真

窗外世界窗內人，入眼隔窗豈分明。

清風本來無跡蹤，月光竹影告動靜。

縱然想悠遊窗外世界，欲得無窗之境，終究受限於『時』、『空』之中的乎！！

凡所有『象』皆是得自『時間』 $T_{begin} \to T_{end}$ 之『取樣』，能得『無窗』耶！！果必『有窗』，可得其『情』乎？？

此理或許太顯然易明，故而 FFT 常常省略其詞？？！！還是有待『無限者』呢！！？？

總需了

窗函數

在信號處理中，窗函數(window function)是一種除在給定區間之外取值均為0的實函數。譬如：在給定區間內為常數而在區間外為0的窗函數被形象地稱為矩形窗。任何函數與窗函數之積仍為窗函數，所以相乘的結果就像透過窗口「看」其他函數一樣。窗函數在頻譜分析、濾波器設計、波束形成、以及音頻數據壓縮（如在Ogg Vorbis音頻格式中）等方面有廣泛的應用。

───

Window function

In signal processing, a window function (also known as an apodization function or tapering function^[1]) is a mathematical function that is zero-valued outside of some chosen interval. For instance, a function that is constant inside the interval and zero elsewhere is called a rectangular window, which describes the shape of its graphical representation. When another function or waveform/data-sequence is multiplied by a window function, the product is also zero-valued outside the interval: all that is left is the part where they overlap, the “view through the window”.

Applications of window functions include spectral analysis, filter design, and beamforming. In typical applications, the window functions used are non-negative, smooth, “bell-shaped” curves.^[2] Rectangle, triangle, and other functions can also be used.

A more general definition of window functions does not require them to be identically zero outside an interval, as long as the product of the window multiplied by its argument is square integrable, and, more specifically, that the function goes sufficiently rapidly toward zero.^[3]

500px-Processing_losses_for_3_window_functions

This figure compares the processing losses of three window functions for sinusoidal inputs, with both minimum and maximum scalloping loss.

───

以及其實務的吧！！！

>>> window = np.hamming(n)
>>> plt.plot(window)
[<matplotlib.lines.Line2D object at 0x30226d0>]
>>> plt.show()

Figure hamming

>>> yw = window * y #Hadamard product
>>> plt.plot(yw)
[<matplotlib.lines.Line2D object at 0x2b0a510>]
>>> plt.show()

Figure sin-hamming

>>> Yw = np.fft.fft(yw)/n
>>> Yw = Yw[range(n/2)]
>>> plt.plot(np.abs(Yw))
[<matplotlib.lines.Line2D object at 0x3034350>]
>>> plt.show()
>>>

Figure sin-w

樹莓派、樹莓派之學習、樹莓派之教育

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】四

2016-05-08 懸鉤子

擊鼓，鼓身靜而鼓面動；操琴，琴弦振然琴柱止。這帶出了『物理現象』之

邊值問題

在微分方程中，邊值問題是一個微分方程和一組稱之為邊界條件的約束條件。邊值問題的解通常是符合約束條件的微分方程的解。

物理學中經常遇到邊值問題，例如波動方程等。許多重要的邊值問題屬於Sturm-Liouville問題。這類問題的分析會和微分算子的本徵函數有關。

在實際應用中，邊值問題應當是適定的（即：存在解，解唯一且解會隨著初始值連續的變化）。許多偏微分方程領域的理論提出是為要證明科學及工程應用的許多邊值問題都是適定問題。

最早研究的邊值問題是狄利克雷問題，是要找出調和函數，也就是拉普拉斯方程的解，後來是用狄利克雷原理找到相關的解。

圖中的區域為微分方程有效的區域，且函數在邊界上的值已知

───

數理之探究關聯上『本徵函數』

Eigenfunction

In mathematics, an eigenfunction of a linear operator D defined on some function space is any non-zero function f in that space that, when acted upon by D, is only multiplied by some scaling factor called an eigenvalue. As an equation, this condition can be written as

Df = \lambda f

for some scalar eigenvalue λ.^[1]^[2]^[3] The solutions to this equation may also be subject to boundary conditions that limit the allowable eigenvalues and eigenfunctions.

An eigenfunction is a type of eigenvector.

……

Applications

Vibrating strings

Let $h (x, t)$ denote the sideways displacement of a stressed elastic chord, such as the vibrating strings of a string instrument, as a function of the position $x$ along the string and of time $t$ . Applying the laws of mechanics to infinitesimal portions of the string, the function $h$ satisfies the partial differential equation

\frac{\partial^2 h}{\partial t^2} = c^2\frac{\partial^2 h}{\partial x^2},

which is called the (one-dimensional) wave equation. Here $c$ is a constant speed that depends on the tension and mass of the string.

This problem is amenable to the method of separation of variables. If we assume that $h (x, t)$ can be written as the product of the form $X (x) T (t)$ , we can form a pair of ordinary differential equations:

\frac{d^2}{dx^2}X=-\frac{\omega^2}{c^2}X, \qquad \frac{d^2}{dt^2}T=-\omega^2 T.

Each of these is an eigenvalue equation with eigenvalues $-\tfrac{\omega^2}{c^2}$ and $- ω 2$ , respectively. For any values of $ω$ and $c$ , the equations are satisfied by the functions

X(x) = \sin\left(\frac{\omega x}{c} + \varphi\right), \qquad T(t) = \sin(\omega t + \psi),

where the phase angles $φ$ and $ψ$ are arbitrary real constants.

If we impose boundary conditions, for example that the ends of the string are fixed at $x = 0$ and $x = L$ , namely $X (0) = X (L) = 0$ , and that $T (0) = 0$ , we constrain the eigenvalues. For these boundary conditions, $sin(φ) = 0$ and $sin(ψ) = 0$ , so the phase angles $φ = ψ = 0$ , and

\sin\left(\frac{\omega L}{c}\right) = 0.

This last boundary condition constrains $ω$ to take a value $ω n = ncπ / L$ , where $n$ is any integer. Thus, the clamped string supports a family of standing waves of the form

h(x,t) = \sin\left(\frac{n\pi x}{L} \right) \sin(\omega_n t).

In the example of a string instrument, the frequency $ω n$ is the frequency of the $n$ ^th harmonic, which is called the $(n - 1)$ ^th overtone.

Standing_wave

The shape of a standing wave in a string fixed at its boundaries is an example of an eigenfunction of a differential operator. The admissible eigenvalues are governed by the length of the string and determine the frequency of oscillation.

───

終於譜出了『傅立葉分析』之美麗花朵︰

駐波形成

在一輛長列『左行』的火車上有一個很長的『水槽』，上有一向右的『行進波』
$u(x, t) = A(x,\ t)\sin (kx - \omega t + \phi)$
，假使向左的火車與向右之水波速度相同，那麼一位站在月台的『觀察者』將如何描述那個『行進波』的呢？

如果觀察水由水龍頭注入水槽的現象，由於水在到達槽底前的流速『較快』，然而到達槽底後水的流速突然的『變慢』，因此會發生『水躍』Hydraulic jump 的現象，此時水之部份動能將轉換為位能，故而在槽底的液面形成『駐波』。這個現象在『河水』的『流速』突然『由快變慢』時也可能發生，因而有人能在『河裡衝浪』，他正站在『駐波』之上！！

那什麼是『駐波』的呢？比方說一個『不動的』stationary 介質中，向左的波 $u_l(k x + \omega t)$ 與向右的波 $u_r(k x - \omega t)$ 疊加後的『合成波』 $u_l +u_r$ ，在『特定』的『邊界條件』下，被『侷限』在一定『空間區域』內無法前進，因此稱為『駐波』。由於駐波不能傳播能量，它的能量將『儲存』在那個空間區域裡。駐波所在區域，『振幅為零』的點稱為『節點』或『波節』Node ，『振幅最大』的點位於兩『節點』之間，通常叫做『腹點』或『波腹』Antinode。

一根長度 $L$ 震盪的弦上，一個向右的簡諧波 $u_r = u_0 \sin(kx - \omega t)$ ，由於弦的兩頭固定，那個波在右端點也只能『反射』回來，形成了 $u_l = u_0 \sin(kx + \omega t)$ ，此時合成波 $u = u_l + u_r$ 是
$u\; = u_0\sin(kx - \omega t) + u_0 \sin(kx + \omega t)$
，可用三角恆等式簡化為
$u = 2 u_0\cos(\omega t)\sin(kx)$
。此時『時間項』與『空間項』分離，形成『駐波』。在 $kx = n \pi$ 時， $\sin(kx) = 0$ ，此處 $n$ 是整數，這就是『節點』；當 $kx = n \pi + \frac{\pi}{2}$ 時 $\parallel \sin(kx) = 1 \parallel$ ，也就是『腹點』。當然波長 $\lambda$ 就得滿足 $\lambda = \frac {L}{n \pi}$ 的邊界條件。

─── 琴弦擇音而振，苟非知音焉得共鳴。───

───

當更能了解那些滿足 $\lambda = \frac {L}{n \pi}$ 『波長』關係的『頻率』構成了那根『弦』的『泛音』。不同『音色』的『弦』正因此『泛音』頻譜不同而出色。或也將知這也是『正交函數族』 Orthogonal functions 的發展以及『傅立葉級數』之歷史濫觴乎︰

Hilbert space interpretation

In the language of Hilbert spaces, the set of functions { $e_n=e^{inx}$ ; n ∈ Z} is an orthonormal basis for the space L²([− $π$ , $π$ ]) of square-integrable functions of [− $π$ , $π$ ]. This space is actually a Hilbert space with an inner product given for any two elements f and g by

\langle f,\, g \rangle \;\stackrel{\mathrm{def}}{=} \; \frac{1}{2\pi}\int_{-\pi}^{\pi} f(x)\overline{g(x)}\,dx.

The basic Fourier series result for Hilbert spaces can be written as

f=\sum_{n=-\infty}^\infty \langle f,e_n \rangle \, e_n.

This corresponds exactly to the complex exponential formulation given above. The version with sines and cosines is also justified with the Hilbert space interpretation. Indeed, the sines and cosines form an orthogonal set:

400px-Fourier_series_integral_identities

Sines and cosines form an orthonormal set, as illustrated above. The integral of sine, cosine and their product is zero (green and red areas are equal, and cancel out) when m, n or the functions are different, and pi only if m and n are equal, and the function used is the same.

\int_{-\pi}^{\pi} \cos(mx)\, \cos(nx)\, dx = \pi \delta_{mn}, \quad m, n \ge 1, \,

\int_{-\pi}^{\pi} \sin(mx)\, \sin(nx)\, dx = \pi \delta_{mn}, \quad m, n \ge 1

(where δ_mn is the Kronecker delta), and

\int_{-\pi}^{\pi} \cos(mx)\, \sin(nx)\, dx = 0;\,

furthermore, the sines and cosines are orthogonal to the constant function 1. An orthonormal basis for L²([− $π$ , $π$ ]) consisting of real functions is formed by the functions 1 and √2 cos(nx), √2 sin(nx) with n = 1, 2,… The density of their span is a consequence of the Stone–Weierstrass theorem, but follows also from the properties of classical kernels like the Fejér kernel.

───

若非如此，樂器將如何和鳴共奏呢？或終可聞箱子天籟之聲的耶？？！！

若問什麼是『頻譜洩漏』

Spectral leakage

The Fourier transform of a function of time, s(t), is a complex-valued function of frequency, S(f), often referred to as a frequency spectrum. Any linear time-invariant operation on s(t) produces a new spectrum of the form H(f)•S(f), which changes the relative magnitudes and/or angles (phase) of the non-zero values of S(f). Any other type of operation creates new frequency components that may be referred to as spectral leakage in the broadest sense. Sampling, for instance, produces leakage, which we call aliases of the original spectral component. For Fourier transform purposes, sampling is modeled as a product between s(t) and a Dirac comb function. The spectrum of a product is the convolution between S(f) and another function, which inevitably creates the new frequency components. But the term ‘leakage’ usually refers to the effect of windowing, which is the product of s(t) with a different kind of function, the window function. Window functions happen to have finite duration, but that is not necessary to create leakage. Multiplication by a time-variant function is sufficient.

Spectral_leakage_from_a_sinusoid_and_rectangular_window

Zoomed view of spectral leakage.

───

？為什麼會有『洩漏』的呢？？奈何讀來宛如天書耶！！或許因為已不記得了史丹佛大學 Brad Osgood 教授開宗明義講︰傅立葉級數假設『週期現象』！！？？

Fourier-1

Fourier-2

───

事實上『取樣數據』 $x_1, x_2, x_3, \cdots , x_k , \cdots x_n$ 之 DFT 只是一種『近似描述』，如果再加上不能滿足『頭尾相銜』、『無縫接續』，有不連續性於其間，必有『效應』的乎？？！！

【看圖說故事】

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> Fs = 150 #取樣頻率
>>> Ts = 1.0/Fs #取樣間隔
>>> t = np.arange(0,1,Ts) #取樣範圍
>>> ff = 5 #訊號頻率
>>> y = np.sin(2 * np.pi * ff * t) #正弦波
>>> plt.subplot(2,1,1)
<matplotlib.axes.AxesSubplot object at 0x2722d90>
>>> plt.plot(t,y,'k-')
[<matplotlib.lines.Line2D object at 0x2387a90>]
>>> plt.xlabel('time')
<matplotlib.text.Text object at 0x2734c50>
>>> plt.ylabel('amplitude')
<matplotlib.text.Text object at 0x2739bd0>
>>> plt.subplot(2,1,2)
<matplotlib.axes.AxesSubplot object at 0x2755a50>
>>> n = len(y)
>>> k = np.arange(n)
>>> T = n/Fs
>>> frq = k/T
>>> freq = frq[range(n/2)]  #單邊頻域
>>> Y = np.fft.fft(y)/n #FFT 且作歸一化
>>> Y = Y[range(n/2)]
>>> plt.plot(freq, abs(Y), 'r-')
[<matplotlib.lines.Line2D object at 0x2755e10>]
>>> plt.xlabel('freq (Hz)')
<matplotlib.text.Text object at 0x2759cd0>
>>> plt.ylabel('|Y(freq)|')
<matplotlib.text.Text object at 0x275cc50>
>>> plt.show()

Figure sin_p

>>> t = np.arange(0,1.1,Ts)
>>> ff = 5
>>> y = np.sin(2 * np.pi * ff * t)
>>> plt.subplot(2,1,1)
<matplotlib.axes.AxesSubplot object at 0x2afcb90>
>>> plt.plot(t,y,'k-')
[<matplotlib.lines.Line2D object at 0x2d67d50>]
>>> plt.xlabel('time')
<matplotlib.text.Text object at 0x2afa290>
>>> plt.ylabel('amplitude')
<matplotlib.text.Text object at 0x2b136d0>
>>> plt.subplot(2,1,2)
<matplotlib.axes.AxesSubplot object at 0x3033950>
>>> n = len(y)
>>> k = np.arange(n)
>>> T = n/Fs
>>> frq = k/T
>>> freq = frq[range(n/2)]
>>> Y = np.fft.fft(y)/n
>>> Y = Y[range(n/2)]
>>> plt.plot(freq, abs(Y), 'r-')
[<matplotlib.lines.Line2D object at 0x3033d10>]
>>> plt.xlabel('freq (Hz)')
<matplotlib.text.Text object at 0x2af9bd0>
>>> plt.ylabel('|Y(freq)|')
<matplotlib.text.Text object at 0x2af3b50>
>>> plt.show()
>>>

Figure sin_l

樹莓派、樹莓派之學習、樹莓派之教育

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】三

2016-05-07 懸鉤子

若是只因『考試！？』而讀書，恐怕會錯過許多學習『樂趣』！！『學問』之道有其『基礎』，『技藝』之術有其『規矩』，所以說本立而道生。若能善用『規矩』，不祇能成『方圓』，又有什麼『圖』不可『製』耶？？

《孟子》‧離婁上
孟子曰：「離婁之明，公輸子之巧，不以規矩，不能成方員：師曠之聰，不以六律，不能正五音；堯舜之道，不以仁政，不能平治天下。今有仁心仁聞而民不被其澤，不可法於後世者，不行先王之道也。

故曰，徒善不足以為政，徒法不能以自行。《詩》云：『不愆不忘，率由舊章。』遵先王之法而過者，未之有也。聖人既竭目力焉，繼之以規矩準繩，以為方員平直，不可勝用也；既竭耳力焉，繼之以六律，正五音，不可勝用也；既竭心思焉，繼之以不忍人之政，而仁覆天下矣。

故曰，為高必因丘陵，為下必因川澤。為政不因先王之道，可謂智乎？是以惟仁者宜在高位。不仁而在高位，是播其惡於眾也。上無道揆也。下無法守也，朝不信道，工不信度，君子犯義，小人犯刑，國之所存者幸也。

故曰，城郭不完，兵甲不多，非國之災也；田野不辟，貨財不聚，非國之害也。上無禮，下無學，賊民興，喪無日矣。《詩》曰：『天之方蹶，無然泄泄。』泄泄，猶沓沓也。事君無義，進退無禮，言則非先王之道者，猶沓沓也。故曰：責難於君謂之恭，陳善閉邪謂之敬，吾君不能謂之賊。

───

將要如何理解『numpy』之 DFT 建置呢︰

There are many ways to define the DFT, varying in the sign of the exponent, normalization, etc. In this implementation, the DFT is defined as

$A_k = \sum_{m=0}^{n-1} a_m \exp\left\{-2\pi i{mk \over n}\right\} \qquad k = 0,\ldots,n-1.$

The DFT is in general defined for complex inputs and outputs, and a single-frequency component at linear frequency $f$ is represented by a complex exponential $a_m = \exp\{2\pi i\,f m\Delta t\}$ , where $\Delta t$ is the sampling interval.

The values in the result follow so-called “standard” order: If A = fft(a, n), then A[0] contains the zero-frequency term (the mean of the signal), which is always purely real for real inputs. Then A[1:n/2] contains the positive-frequency terms, and A[n/2+1:] contains the negative-frequency terms, in order of decreasingly negative frequency. For an even number of input points, A[n/2] represents both positive and negative Nyquist frequency, and is also purely real for real input. For an odd number of input points, A[(n-1)/2] contains the largest positive frequency, while A[(n+1)/2] contains the largest negative frequency. The routine np.fft.fftfreq(n) returns an array giving the frequencies of corresponding elements in the output. The routine np.fft.fftshift(A) shifts transforms and their frequencies to put the zero-frequency components in the middle, and np.fft.ifftshift(A) undoes that shift.

When the input a is a time-domain signal and A = fft(a), np.abs(A) is its amplitude spectrum and np.abs(A)**2 is its power spectrum. The phase spectrum is obtained by np.angle(A).

The inverse DFT is defined as

$a_m = \frac{1}{n}\sum_{k=0}^{n-1}A_k\exp\left\{2\pi i{mk\over n}\right\} \qquad m = 0,\ldots,n-1.$

It differs from the forward transform by the sign of the exponential argument and the default normalization by $1/n$ .

───

即使輔之以維基百科詞條說明︰

離散傅立葉變換

離散傅立葉變換（Discrete Fourier Transform，縮寫為DFT），是傅立葉變換在時域和頻域上都呈離散的形式，將信號的時域採樣變換為其DTFT的頻域採樣。在形式上，變換兩端（時域和頻域上）的序列是有限長的，而實際上這兩組序列都應當被認為是離散周期信號的主值序列。即使對有限長的離散信號作DFT，也應當將其看作其周期延拓的變換。在實際應用中通常採用快速傅立葉變換計算DFT。

定義

對於N點序列 $\left\{x[n]\right\}_{0\le n <N}$ ，它的離散傅立葉變換（DFT）為

\hat{x}[k]=\sum_{n=0}^{N-1} e^{-i\frac{2\pi}{N}nk}x[n] \qquad k = 0,1,\ldots,N-1.

其中 $e$ 是自然對數的底數， $i$ 是虛數單位。通常以符號 $\mathcal{F}$ 表示這一變換，即

\hat{x}=\mathcal{F}x

離散傅立葉變換的逆變換（IDFT）為：

x\left[n\right]={1 \over N}\sum_{k=0}^{N-1} e^{ i\frac{2\pi}{N}nk}\hat{x}[k] \qquad n = 0,1,\ldots,N-1.

可以記為：

x=\mathcal{F}^{-1}\hat{x}

實際上，DFT和IDFT變換式中和式前面的歸一化係數並不重要。在上面的定義中，DFT和IDFT前的係數分別為1和1/N。有時會將這兩個係數都改成 $1/\sqrt{N}$ 。

從連續到離散

連續時間信號x(t)以及對應的連續傅立葉變換 $\hat{x}(\omega)$ 都是連續函數。由於數字系統只能處理有限長的離散信號，因此必須將 $x$ 和 $\hat{x}$ 都離散化，並且建立對應的傅立葉變換。

假設x（t）時限於[0, L]，再通過時域採樣將 $x(t)$ 離散化，就可以得到有限長離散信號，記為 $x_{discrete}(t)$ 。設採樣周期為T，則時域採樣點數N=L/T。

x_{discrete}(t) = x(t)\sum_{n=0}^{N-1}\delta(t-nT)=\sum_{n=0}^{N-1}x(nT)\delta(t-nT)

它的傅立葉變換為

\hat{x}_{discrete}(\omega) = \sum_{n=0}^{N-1}x(nT)\mathcal{F}\delta(t-nT) = \sum_{n=0}^{N-1}x(nT)e^{-i n\omega T}

這就是x(t)在時域採樣後的連續傅立葉變換，也就是離散時間傅立葉變換，它在頻域依然是連續的。

下面將頻域信號轉化為有限長離散信號。與對時域信號的處理類似，假設頻域信號是帶限的，再經過離散化，即可得到有限長離散信號。依據採樣定理，時域採樣若要能完全重建原信號，頻域信號 $\hat{x}(\omega)$ 應當帶限於（0,1/（2*T））。由於時域信號時限於[0, L]，由採樣定理以及時頻對偶的關係，頻域的採樣間隔應為1/L。故，頻域採樣點數為：

\frac{1/T}{1/L}=N

即頻域採樣的點數和時域採樣同為N，頻域採樣點為 $\{\omega_k={2\pi}k/NT\}_{0 \le k < N}$ 在DTFT頻域上採樣：

\hat{x}[k] = \hat{x}_{discrete}(\omega_k) = \frac1{T} \sum_{n=0}^{N-1}x[nT]e^{-i\frac{2\pi}{N} nk}

令T=1，將其歸一化，就得到前面定義的離散傅立葉變換。因此，DFT就是先將信號在時域離散化，求其連續傅立葉變換後，再在頻域離散化的結果。

───

雖有助益，或許難以通熟。

何不好好了解『複數』是什麼？？

那麽要怎樣理解『複數』 $z = x + i \ y$ 的呢？如果說『複數』起源於『方程式』的『求解』，比方說 $x^2 + 1 = 0, \ x = \pm i$ ，這定義了『 $i = \sqrt{-1}$ 』，但是它的『意義』依然晦澀。即使說從『複數平面』的每一個『點』 $(x, y)$ 都對應著一個『複數』 $z = x + i \ y$ 可能還是不清楚『 $i$ 』的意思到底是什麼？假使再從『複數』的『加法上看』︰

假使 $z_1 = x_1 + i \ y_1$ 、 $z_2 = x_2 + i \ y_2$ ，

那麼 $z_1 + z_2 = (x_1 + x_2) + i \ (y_1 + y_2)$ 。

這是一種類似『向量』的加法，是否『 $i$ 』的意義就藏在其中的呢？

一九九八年美國新罕布希爾大學 University of New Hampshire 的
Paul J. Nahin 教授寫了一本『An Imaginary Tale: the Story of the Square Root of −1』的書，指出韋塞爾當初所講的『幾何意義』就是︰

$i = \sqrt{-1} = 1 \ \angle 90^{\circ}$

，也就是說『i』就是『逆時鐘旋轉九十度』的『運算子』！

假使從複數的『極座標』表示法來看複數的『乘法』︰

假使 $z_1 = r \cdot e^{i \ \theta}, \ z_2 = \alpha \cdot e^{i \ \beta}$ ，那麼 $z_1 \cdot z_2 = \alpha \cdot r \cdot e^{i \ (\theta +\beta)}$

。就可以解釋成 Z1 『向量』被『逆時鐘旋轉』了『β』角度，它的『長度』被『縮放』了『α』倍！！

複數果真不是簡單的『數』啊！也難怪它是『完備的』的喔！！

電子和工程領域中，常常會使用到『正弦』 Sin 信號，一般可以使用『相量』 Phasor 來作簡化分析。『相量』是一個『複數』，也是一種『向量』，通常使用『極座標』表示，舉例來說一個『振幅』是 $A$ ，『角頻率』是 $\omega$ ，初始『相位角』是 $\theta$ 的『正弦信號』可以表示為 $A \ e^{j \ (\omega t + \theta)}$ ，這裡的『j』就是『複數的 i』。為什麼又要改用 $j = \sqrt{-1}$ 的呢？這是因為再『電子學』和『電路學』領域中 $i$ 通常代表著『電流』， $v$ 通常代表了『電壓』，因此為了避免『混淆』起見，所以才會『更名用 j』。

尤拉公式 Euler’s formula，是複數分析中的公式，它將三角函數與複數指數函數相關聯，對任意實數 $x$ ，都有

$e^{j x} = \cos x + j \sin x$

，它的重要性是不言而喻的啊！！

─── 摘自《【Sonic π】電聲學補充《二》》

踏踏實實的上一堂課！！

若想更深入的了解，何不上一堂 Stanford 大學的公開課︰

Stanford Engineering Everywhere EE261 – The Fourier Transform and its Applications

author: Brad G. Osgood, Computer Science Department, Stanford University
released under terms of: Creative Commons Attribution Non-Commercial (CC-BY-NC)

The goals for the course are to gain a facility with using the Fourier transform, both specific techniques and general principles, and learning to recognize when, why, and how it is used. Together with a great variety, the subject also has a great coherence, and the hope is students come to appreciate both.

……

同時認真讀讀 Brad Osgood 教授之課堂筆記耶。

Lecture Notes for

EE 261

The Fourier Transform and its Applications

Prof. Brad Osgood
Electrical Engineering Department
Stanford University

FT-1

FT-2

─── 摘自《勇闖新世界︰ W!o《卡夫卡村》變形祭︰品味科學‧教具教材‧【專題】 PD‧箱子世界‧傅立葉》

『真積力，久則入』，一朝自然水到磲成的乎！！？？

樹莓派、樹莓派之學習、樹莓派之教育

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】二

2016-05-06 懸鉤子

什麼是事物的『特徵』呢？為什麼它的『提取方法』很重要？維基百科詞條這麼說︰

Feature extraction

In machine learning, pattern recognition and in image processing, feature extraction starts from an initial set of measured data and builds derived values (features) intended to be informative and non-redundant, facilitating the subsequent learning and generalization steps, and in some cases leading to better human interpretations. Feature extraction is related to dimensionality reduction.

When the input data to an algorithm is too large to be processed and it is suspected to be redundant (e.g. the same measurement in both feet and meters, or the repetitiveness of images presented as pixels), then it can be transformed into a reduced set of features (also named a features vector). This process is called feature selection. The selected features are expected to contain the relevant information from the input data, so that the desired task can be performed by using this reduced representation instead of the complete initial data.

───

假使考慮如何『定義』事物耶？也許『特徵』就是『界定性徵』，可以用來『區分』相異的東西！所以人們自然懂得『汪星人』不同於『喵星人』的也！！

於是乎好奇那『聲音』本有『調子』，可以用

Cepstrum

A cepstrum (/ˈkɛpstrəmˈˌˈsɛpstrəmˈ/) is the result of taking the Inverse Fourier transform (IFT) of the logarithm of the estimated spectrum of a signal. It may be pronounced in the two ways given, the second having the advantage of avoiding confusion with ‘kepstrum’ which also exists (see below). There is a complex cepstrum, a real cepstrum, a power cepstrum, and a phase cepstrum. The power cepstrum in particular finds applications in the analysis of human speech.

The name “cepstrum” was derived by reversing the first four letters of “spectrum”. Operations on cepstra are labelled quefrency analysis (aka quefrency alanysis^[1]), liftering, or cepstral analysis.

Cepstrum_signal_analysis

Steps in forming cepstrum from time history

───

來探討。那麼『圖象』可有『調子』乎？！能否依樣畫葫蘆來研究的呢！？不管『笨鳥先飛』、『菜鳥忘飛』、『老鳥已飛』……… 科技史裡滿載『傻問題』之『大成就』矣！！？？何不就效法一下嘛？？！！

【還是用五】

>>> img = training_data[0][0].reshape(28,28)
>>> f_img = network.np.fft.rfft2(img)
>>> logp_img = 2*network.np.log(network.np.abs(f_img))
>>> plt.imshow(logp_img)
<matplotlib.image.AxesImage object at 0x51af290>
>>> plt.show()

Figure 5p

>>> ilogpf_img = network.np.fft.irfft2(logp_img)
>>> cf_img = network.np.abs(ilogpf_img)**2
>>> plt.imshow(cf_img)
<matplotlib.image.AxesImage object at 0x51d1050>
>>> plt.show()

Figure 5c

【依舊選零】

>>> img1 = training_data[1][0].reshape(28,28)
>>> f_img1 = network.np.fft.rfft2(img1)
>>> logp_img1 = 2*network.np.log(network.np.abs(f_img1))
>>> plt.imshow(logp_img1)
<matplotlib.image.AxesImage object at 0x5091e50>
>>> plt.show()

Figure 0p

>>> ilogpf_img1 = network.np.fft.irfft2(logp_img1)
>>> cf_img1 = network.np.abs(ilogpf_img1)**2
>>> plt.imshow(cf_img1)
<matplotlib.image.AxesImage object at 0x51b9b50>
>>> plt.show()
>>>

Figure 0c

【參考資料】

Discrete Fourier Transform (`numpy.fft`)

Standard FFTs

`fft`(a[, n, axis, norm])	Compute the one-dimensional discrete Fourier Transform.
`ifft`(a[, n, axis, norm])	Compute the one-dimensional inverse discrete Fourier Transform.
`fft2`(a[, s, axes, norm])	Compute the 2-dimensional discrete Fourier Transform This function computes the n-dimensional discrete Fourier Transform over any axes in an M-dimensional array by means of the Fast Fourier Transform (FFT).
`ifft2`(a[, s, axes, norm])	Compute the 2-dimensional inverse discrete Fourier Transform.
`fftn`(a[, s, axes, norm])	Compute the N-dimensional discrete Fourier Transform.
`ifftn`(a[, s, axes, norm])	Compute the N-dimensional inverse discrete Fourier Transform.

Real FFTs

`rfft`(a[, n, axis, norm])	Compute the one-dimensional discrete Fourier Transform for real input.
`irfft`(a[, n, axis, norm])	Compute the inverse of the n-point DFT for real input.
`rfft2`(a[, s, axes, norm])	Compute the 2-dimensional FFT of a real array.
`irfft2`(a[, s, axes, norm])	Compute the 2-dimensional inverse FFT of a real array.
`rfftn`(a[, s, axes, norm])	Compute the N-dimensional discrete Fourier Transform for real input.
`irfftn`(a[, s, axes, norm])	Compute the inverse of the N-dimensional FFT of real input.

Hermitian FFTs

`hfft`(a[, n, axis, norm])	Compute the FFT of a signal which has Hermitian symmetry (real spectrum).
`ihfft`(a[, n, axis, norm])	Compute the inverse FFT of a signal which has Hermitian symmetry.

Helper routines

`fftfreq`(n[, d])	Return the Discrete Fourier Transform sample frequencies.
`rfftfreq`(n[, d])	Return the Discrete Fourier Transform sample frequencies (for usage with rfft, irfft).
`fftshift`(x[, axes])	Shift the zero-frequency component to the center of the spectrum.
`ifftshift`(x[, axes])	The inverse of `fftshift`.

Background information

Fourier analysis is fundamentally a method for expressing a function as a sum of periodic components, and for recovering the function from those components. When both the function and its Fourier transform are replaced with discretized counterparts, it is called the discrete Fourier transform (DFT). The DFT has become a mainstay of numerical computing in part because of a very fast algorithm for computing it, called the Fast Fourier Transform (FFT), which was known to Gauss (1805) and was brought to light in its current form by Cooley and Tukey [CT]. Press et al. [NR] provide an accessible introduction to Fourier analysis and its applications.

Because the discrete Fourier transform separates its input into components that contribute at discrete frequencies, it has a great number of applications in digital signal processing, e.g., for filtering, and in this context the discretized input to the transform is customarily referred to as a signal, which exists in the time domain. The output is called a spectrum or transform and exists in the frequency domain.

Implementation details

There are many ways to define the DFT, varying in the sign of the exponent, normalization, etc. In this implementation, the DFT is defined as

$A_k = \sum_{m=0}^{n-1} a_m \exp\left\{-2\pi i{mk \over n}\right\} \qquad k = 0,\ldots,n-1.$

When the input a is a time-domain signal and A = fft(a), np.abs(A) is its amplitude spectrum and np.abs(A)**2 is its power spectrum. The phase spectrum is obtained by np.angle(A).

The inverse DFT is defined as

$a_m = \frac{1}{n}\sum_{k=0}^{n-1}A_k\exp\left\{2\pi i{mk\over n}\right\} \qquad m = 0,\ldots,n-1.$

It differs from the forward transform by the sign of the exponential argument and the default normalization by $1/n$ .

───

竟然會看起來很像？？似乎又有點不一樣！！到底該說是『行』還是『不行』的呀？？？

FreeSandal

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】六

短時距傅立葉變換

與傅立葉轉換在概念上的區別

定義

數學定義

窗函數

頻譜(Spectrogram)

Short-time Fourier transform

Resolution issues

Fourier transform on Euclidean space

Uncertainty principle

scipy.signal.spectrogram

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】五

窗函數

Window function

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】四

邊值問題

Eigenfunction

Applications

Vibrating strings

Hilbert space interpretation

Spectral leakage

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】三

離散傅立葉變換

定義

從連續到離散

Stanford Engineering Everywhere EE261 – The Fourier Transform and its Applications

W!o+ 的《小伶鼬工坊演義》︰神經網絡【FFT】二

Feature extraction

Cepstrum

Discrete Fourier Transform (`numpy.fft`)

Standard FFTs

Real FFTs

Hermitian FFTs

Helper routines

Background information

Implementation details

輕。鬆。學。部落客

2024 年 5 月
日	一	二	三	四	五	六
« 4 月
			1	2	3	4
5	6	7	8	9	10	11
12	13	14	15	16	17	18
19	20	21	22	23	24	25
26	27	28	29	30	31

與傅立葉轉換在概念上的區別

定義

數學定義

窗函數

頻譜(Spectrogram)

Resolution issues

Uncertainty principle

Applications

Vibrating strings

定義

從連續到離散

Discrete Fourier Transform (numpy.fft)

Standard FFTs

Real FFTs

Hermitian FFTs

Helper routines

Background information

Implementation details

輕。鬆。學。部落客

Discrete Fourier Transform (`numpy.fft`)