関数morlet_spectrumマニュアル

(The documentation of function morlet_spectrum)

Last Update: 2025/8/7

◆機能・用途(Purpose)

Morlet wavelet \[\begin{equation} \psi(t) =\pi^{-1/4}\left[\exp\left(-2\pi if_{ref} t\right) -\exp\left(-2\pi^2 f_{ref}^2\right)\right] \exp\left(-\frac{t^2}{2}\right) \label{eq.psi} \end{equation}\] のフーリエ変換を計算する。求めるのは指定した1つの周波数での値とし、以下の(\ref{eq.Psi.final})式を用いて計算する。
Compute the Fourier transformation of a Morlet wavelet (Eq. \ref{eq.psi}) for a single, specified frequency using Eq. \ref{eq.Psi.final}.

◆形式(Format)

#include <sequence/wavelet.h>
inline double morlet_spectrum(const double f,const double f_ref)

◆引数(Arguments)

f	フーリエ変換の値を求めたい周波数(Hz)。以下の計算式における\(f\)。 The frequency (Hz) at which the Fourier transformation is needed, which is the value \(f\) in the formula below.
f_ref	(\ref{eq.psi})式における\(f_{ref}\)の値。 The value of \(f_{ref}) in Eq. (\ref{eq.psi}).

◆戻り値(Return value)

指定した周波数\(f\)におけるフーリエ変換\(\Psi(f)\)の値。
The value of Fourier transformation \(\Psi(f)\) at the specified frequency \(f\).

◆使用例(Example)

double Psi=morlet_spectrum(5.1,3.5);

◆計算式(Formula)

sequence/fft.hに合わせてフーリエ変換を \[\begin{equation} \Psi(f)\equiv\int_{-\infty}^{\infty}\psi(t)\exp(2\pi ift)dt \label{eq.psi.fourier} \end{equation}\] と定義する。(\ref{eq.psi.fourier})に(\ref{eq.psi})を代入すると \[\begin{eqnarray} \Psi(f) &=& \int_{-\infty}^{\infty} \pi^{-1/4}\left[\exp\left(-2\pi if_{ref} t\right) -\exp\left(-2\pi^2 f_{ref}^2\right)\right] \exp\left(-\frac{t^2}{2}\right)\exp(2\pi ift)dt \nonumber \\ &=& \pi^{-1/4}\int_{-\infty}^{\infty} \exp\left[2\pi i\left(f-f_{ref}\right)t-\frac{t^2}{2}\right]dt -\pi^{-1/4}\exp\left(-2\pi^2 f_{ref}^2\right) \nonumber \\ & & \int_{-\infty}^{\infty}\exp\left[2\pi ift-\frac{t^2}{2}\right]dt \nonumber \\ &=& \pi^{-1/4}\int_{-\infty}^{\infty} \exp\left[-\frac{t^2-4\pi i\left(f-f_{ref}\right)t}{2}\right]dt -\pi^{-1/4}\exp\left(-2\pi^2 f_{ref}^2\right) \nonumber \\ & & \int_{-\infty}^{\infty}\exp\left[-\frac{t^2-4\pi ift}{2}\right]dt \nonumber \\ &=& \pi^{-1/4}\int_{-\infty}^{\infty} \exp\left[-\frac{t^2-4\pi i\left(f-f_{ref}\right)t +4\pi^2 i^2\left(f-f_{ref}\right)^2 -4\pi^2 i^2\left(f-f_{ref}\right)^2}{2}\right]dt \nonumber \\ & & -\pi^{-1/4}\exp\left(-2\pi^2 f_{ref}^2\right) \int_{-\infty}^{\infty} \exp\left[-\frac{t^2-4\pi ift+4\pi^2 i^2f^2-4\pi^2 i^2f^2}{2}\right]dt \nonumber \\ &=& \pi^{-1/4}\int_{-\infty}^{\infty} \exp\left[-\frac{\left\{t-2\pi i\left(f-f_{ref}\right)\right\}^2}{2} -2\pi^2 \left(f-f_{ref}\right)^2\right]dt \nonumber \\ & & -\pi^{-1/4}\exp\left(-2\pi^2 f_{ref}^2\right) \int_{-\infty}^{\infty} \exp\left[-\frac{(t-2\pi if)^2}{2}-2\pi^2 f^2\right]dt \nonumber \\ &=& \pi^{-1/4}\exp\left[-2\pi^2 \left(f-f_{ref}\right)^2\right] \int_{-\infty}^{\infty} \exp\left[-\frac{\left\{t-2\pi i\left(f-f_{ref}\right)\right\}^2}{2} \right]dt \nonumber \\ & & -\pi^{-1/4}\exp\left[-2\pi^2 \left(f^2+f_{ref}^2\right)\right] \int_{-\infty}^{\infty} \exp\left[-\frac{(t-2\pi if)^2}{2}\right]dt \label{eq.Psi.derivation1} \end{eqnarray}\] となる。ここで複素平面上の任意の2地点\(P\), \(Q\)を結ぶ直線に沿った複素積分 \[\begin{equation} I_{P\rightarrow Q}(f)\equiv\int_{P\rightarrow Q} \exp\left[-\frac{(z-2\pi if)^2}{2}\right]dz \label{eq.integral.complex.general} \end{equation}\] を考える。被積分関数は複素平面全体にわたって正則であるので、 \(T\)を正の実数として複素平面上の4地点 \(A(T)\equiv -T\), \(B(T)\equiv T\), \(C(T,f)\equiv T+2\pi if\), \(D(T,f)\equiv -T+2\pi if\) を結んだ周回積分は \[\begin{equation} I_{A(T)\rightarrow B(T)}(f) +I_{B(T)\rightarrow C(T,f)}(f) +I_{C(T,f)\rightarrow D(T,f)}(f) +I_{D(T,f)\rightarrow A(T,f)}(f) =0 \label{eq.round_integral} \end{equation}\] となる。ここで\(T\rightarrow\infty\)とした極限を考えると

\(A(T)\)から\(B(T)\)までの積分(\(z=t\), \(dz=dt\)): \[\begin{eqnarray} I_{A(T)\rightarrow B(T)}(f) &=& \int_{-T}^T \exp\left[-\frac{(t-2\pi if)^2}{2}\right]dt \nonumber \\ &\rightarrow& \int_{-\infty}^{\infty} \exp\left[-\frac{(t-2\pi if)^2}{2}\right]dt \label{eq.integral.AB} \end{eqnarray}\]
\(B(T)\)から\(C(T,f)\)までの積分(\(z=T+is\), \(dz=ids\)): \[\begin{eqnarray*} I_{B(T)\rightarrow C(T,f)}(f) &=& \int_{0}^{2\pi f}\exp\left[-\frac{(T+is-2\pi if)^2}{2}\right]ids \nonumber \\ &=& i\int_{0}^{2\pi f} \exp\left[-\frac{T^2+2iT(s-2\pi f)+i^2(s-2\pi f)^2}{2}\right]ds \nonumber \\ &=& i\exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \exp\left[-iT(s-2\pi f)\right]\exp\left[\frac{(s-2\pi f)^2}{2}\right]ds \end{eqnarray*}\] \[\begin{eqnarray} \left|I_{B(T)\rightarrow C(T,f)}(f)\right| &=& \exp\left(-\frac{T^2}{2}\right)\left|\int_{0}^{2\pi f} \exp\left[-iT(s-2\pi f)\right]\exp\left[\frac{(s-2\pi f)^2}{2}\right]ds \right| \nonumber \\ &\leq& \exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \left|\exp\left[-iT(s-2\pi f)\right] \exp\left[\frac{(s-2\pi f)^2}{2}\right]\right|ds \nonumber \\ &=& \exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \exp\left[\frac{(s-2\pi f)^2}{2}\right]ds \nonumber \\ &\leq& \exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \exp\left[\frac{(-2\pi f)^2}{2}\right]ds \nonumber \\ &=& \exp\left(-\frac{T^2}{2}\right) \exp\left[\frac{(-2\pi f)^2}{2}\right]\cdot 2\pi f \nonumber \\ &\rightarrow& 0 \label{eq.integral.BC} \end{eqnarray}\]
\(C(T,f)\)から\(D(T,f)\)までの積分(\(z=t+2\pi if\), \(dz=dt\)): \[\begin{eqnarray} I_{C(T,f)\rightarrow D(T,f)}(f) &=& \int_{T}^{-T} \exp\left[-\frac{(t+2\pi if-2\pi if)^2}{2}\right]dt \nonumber \\ &=& -\int_{-T}^{T} \exp\left[-\frac{t^2}{2}\right]dt \nonumber \\ &\rightarrow& -\int_{-\infty}^{\infty}\exp\left[-\frac{t^2}{2}\right]dt \nonumber \\ &=& -\sqrt{2\pi} \label{eq.integral.CD} \end{eqnarray}\]
\(D(T,f)\)から\(A(T)\)までの積分(\(z=-T+is\), \(dz=ids\)): \[\begin{eqnarray*} I_{D(T,f)\rightarrow A(T)}(f) &=& \int_{2\pi f}^{0}\exp\left[-\frac{(-T+is-2\pi if)^2}{2}\right]ids \nonumber \\ &=& -i\int_{0}^{2\pi f} \exp\left[-\frac{T^2-2iT(s-2\pi f)+i^2(s-2\pi f)^2}{2}\right]ds \nonumber \\ &=& -i\exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \exp\left[iT(s-2\pi f)\right]\exp\left[\frac{(s-2\pi f)^2}{2}\right]ds \end{eqnarray*}\] \[\begin{eqnarray} \left|I_{D(T,f)\rightarrow A(T)}(f)\right| &=& \exp\left(-\frac{T^2}{2}\right)\left|\int_{0}^{2\pi f} \exp\left[iT(s-2\pi f)\right]\exp\left[\frac{(s-2\pi f)^2}{2}\right]ds \right| \nonumber \\ &\leq& \exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \left|\exp\left[iT(s-2\pi f)\right] \exp\left[\frac{(s-2\pi f)^2}{2}\right]\right|ds \nonumber \\ &=& \exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \exp\left[\frac{(s-2\pi f)^2}{2}\right]ds \nonumber \\ &\leq& \exp\left(-\frac{T^2}{2}\right)\int_{0}^{2\pi f} \exp\left[\frac{(-2\pi f)^2}{2}\right]ds \nonumber \\ &=& \exp\left(-\frac{T^2}{2}\right) \exp\left[\frac{(-2\pi f)^2}{2}\right]\cdot 2\pi f \nonumber \\ &\rightarrow& 0 \label{eq.integral.DA} \end{eqnarray}\]

となるので、(\ref{eq.round_integral})式で\(T\rightarrow\infty\)とした極限は \[\begin{equation} \int_{-\infty}^{\infty}\exp\left[-\frac{(t-2\pi if)^2}{2}\right]dt -\sqrt{2\pi}=0 \label{eq.round_integral.infT} \end{equation}\] となり、ここから \[\begin{equation} \int_{-\infty}^{\infty}\exp\left[-\frac{(t-2\pi if)^2}{2}\right]dt =\sqrt{2\pi} \label{eq.integral1} \end{equation}\] が得られる。同様にして \[\begin{equation} \int_{-\infty}^{\infty} \exp\left[-\frac{\left\{t-2\pi i\left(f-f_{ref}\right)\right\}^2}{2} \right]dt =\sqrt{2\pi} \label{eq.integral2} \end{equation}\] も得られ、これらを(\ref{eq.Psi.derivation1})に代入すると \[\begin{eqnarray} \Psi(f) &=& \pi^{-1/4}\exp\left[-2\pi^2 \left(f-f_{ref}\right)^2\right] \sqrt{2\pi} -\pi^{-1/4}\exp\left[-2\pi^2 \left(f^2+f_{ref}^2\right)\right] \sqrt{2\pi} \nonumber \\ &=& \sqrt{2}\pi^{1/4}\exp\left[-2\pi^2 \left(f-f_{ref}\right)^2\right] -\sqrt{2}\pi^{1/4}\exp\left[-2\pi^2 \left(f^2+f_{ref}^2\right)\right] \label{eq.Psi.final} \end{eqnarray}\] が得られる。
Following sequence/fft.h, the Fourier transformation is defined by Eq. (\ref{eq.psi.fourier}). Inserting Eq. (\ref{eq.psi}) into (\ref{eq.psi.fourier}) results in (\ref{eq.Psi.derivation1}). Let us consider a complex integral (\ref{eq.integral.complex.general}) along a straight line connecting two arbitrary points \(P\) and \(Q\) on a complex plane. Since the integrand is regular over the entire complex plane, a round integral connecting 4 points \(A(T)\equiv -T\), \(B(T)\equiv T\), \(C(T,f)\equiv T+2\pi if\), and \(D(T,f)\equiv -T+2\pi if\) on a complex plane, where \(T\) is a positive real number, satisfies (\ref{eq.round_integral}). For a limit of \(T\rightarrow\infty\), each term in Eq. (\ref{eq.round_integral}) can be calculated as follows:

the integral from \(A(T)\) to \(B(T)\), where \(z=t\) and \(dz=dt\), is calculated as (\ref{eq.integral.AB});
the integral from \(B(T)\) to \(C(T,f)\), where \(z=T+is\) and \(dz=ids\), is calculated as (\ref{eq.integral.BC});
the integral from \(C(T,f)\) to \(D(T,f)\), where \(z=t+2\pi if\) and \(dz=dt\), is calculated as (\ref{eq.integral.CD}); and
the integral from \(D(T,f)\) to \(A(T)\), where \(z=-T+is\) and \(dz=ids\), is calculated as (\ref{eq.integral.DA}).

Therefore Eq. (\ref{eq.round_integral}) reduces to (\ref{eq.round_integral.infT}) at the limit of \(T\rightarrow\infty\), giving (\ref{eq.integral1}). In the same way, Eq. (\ref{eq.integral2}) is obtained. Inserting these results into (\ref{eq.Psi.derivation1}) results in (\ref{eq.Psi.final}).

◆検証(Validation)

(\ref{eq.psi})式に基づいてMorlet waveletを生成して数値的にフーリエ変換し、本関数で計算した値と比較した。
I compared the numerical Fourier transformation of a Morlet wavelet created by Eq. (\ref{eq.psi}) with a spectrum created by this function.

Morlet waveletの生成にあたってはデータサンプル数\(N=16384\)、時間刻み\(\Delta t=0.01\)、先頭時刻\(t_0=-N\Delta t/2\)とした。時刻0付近に最も強いエネルギーを持つ波形なので先頭時刻ではなく中心時刻が0になるように時刻範囲を設定することがポイントである。 Morlet waveletは複素関数なのでstruct imsequence2型構造体に代入し、これを関数fft2によってスペクトルに変換した。
To generate the Morlet wavelet, I used the number of data samples \(N=16384\), a time stepping \(\Delta t=0.01\), and a start time \(t_0=-N\Delta t/2\). The point here is to set the time range not starting at time 0 but centered on time 0, because the wavelet has the strongest energy at the time 0. The Morlet wavelet is inserted into a struct imsequence2-type structure since it is a complex function. Then the wavelet is transformed to a spectrum by function fft2.

このようにして生成したスペクトルは実部に比べて虚部が十分に小さく、実関数と見なせるものであった。またその実部はこの関数を用いて計算したスペクトルとグラフ上で差異が認められないレベルで一致した。なお、この検証は\(f_{ref}=1\)および\(f_{ref}=3.3\)の2つのケースで試した。
The spectrum created in this way had the imaginary components smaller enough than the real ones to regard the spectrum as a real function. The real parts showed no difference at a visible scale in a graph from the spectrum calculated with this function. This examination was conducted for two cases, namely, \(f_{ref}=1\) and \(f_{ref}=3.3\).