Algorithm of Wavelet Transform

There are many kind of wavelet transform, and I think the names are quite confusing.

For instance, there are continuous and discrete wavelet transforms, in which, the “continuous” and “discrete” are for the wavelet parameters, not for the “data” itself. Therefore, for discrete data, there are “continuous” and “discrete” wavelet transforms, and for function, there are also “continuous” and “discrete” wavelet transforms.

In here, we will focus on discrete wavelet transform for function first. This discrete wavelet transform is also called as wavelet series, which express a compact support function into series of wavelet.

For simplicity, we also focus on orthonormal wavelet.

As the wavelet span the entire space, any compact function can be expressed as

$f(t) = \sum_{j,k} \left \psi_{j,k}(t)$

$\psi_{j,k}(t) = 2^{j/2} \psi(2^j t - k)$

where $j, k$ are integer.

Now, we move to discrete data discrete wavelet transform. The data is discrete, we can imagine only $t_n = t_0 + n \Delta$ points are known with finite $n$.

$f_n = f(t_n) = \sum_{j,k} \left \psi_{j,k}(t_n)$

the integration becomes a finite sum.

Without loss of generality, we can set $t_0 = 0, \Delta = 1$, and then the time axis becomes an integer number axis. We found that $j \leq 0$ as the wavelet can only be expand, not shrink. Because there are finite number of data point, i.e. $n < \infty$, $-Log_2(n) < j \leq 0$.

However, this double summation for each $f_n$ is very time consuming. There is a Fast Discrete Wavelet Transform. Before we continuous, we must study the wavelet.

From the last post, we know that the scaling function that generate a MRA must be:

$\phi(t) = \sum_{k} g_0(k) \phi(2t-k)$

$\left<\phi(t-k) | \phi(t-k') \right> = \delta_{kk'}$

, where $k$ are integer. The set of shifted scaling function span a space $V_0$. For the wavelet,

$\psi(t) = \sum_{k} g_1(k) \psi(2t-k)$

$\left<\psi(t-k) | \psi(t-k') \right> = \delta_{kk'}$

The set of shifted wavelet span a space $W_0$, so that $W_0 \perp V_0$, so that

$\left<\phi(t-k)|\psi(t-k') \right> = 0$

Since the wavelet is generated from the scaling function, we expect the coefficient of $g_0(k)$ and $g_1(k)$ are related. In fact, the relationship for orthonormal scaling function and wavelet is

$g_1(k) = (-1)^k g_0(1-k)$

I will continuous later….. very tired….

Coulomb Energy in Nuclear Unit

Coulomb energy between 2 charge $q_1, q_2$ in SI unit is

$\displaystyle U_c = \frac{1}{4\pi\epsilon_0} \frac{q_1 q_2}{r} [J]$

$\displaystyle U_c = \frac{e^2}{4\pi\epsilon_0} \frac{Z_1 Z_2}{r [m]} [J \cdot m] = \frac{Z_1 Z_2}{r [m]} (\frac{e^2}{4\pi \epsilon_0}) [J \cdot m] = \frac{Z_1 Z_2}{r [m]} 2.30708 \times 10^{-28} [J \cdot m]$

We need to convert the SI unit into nuclear unit:

$1 [J] = \frac{10^{-6}}{e} [MeV]$

$1 [m] = 10^{15} [fm]$

$\displaystyle U_c = \frac{Z_1 Z_2}{r [m]} 1.43996 \times 10^{-09} [MeV \cdot m]$

$\displaystyle U_c = \frac{Z_1 Z_2}{r [fm]} 1.43996 [MeV \cdot fm]$

Therefore, a simple expression

$\displaystyle U_c = \frac{e^2}{r} Z_1 Z_2$

where $e^2 = 1.44 [MeV\cdot fm]$

Other useful quantities are:

• $\hbar c = 197.327 [Mev\cdots fm]$
• $e^2/\hbar c = 1/137.036$
• $\hbar = 6.58212 [MeV\cdot s]$
• $c = 2.99792458\times 10^{23} [fm/s]$

Qt with FFTW

To use Qt with FFTW (Fast Fourier Transfrom in the West) in windows.

http://www.fftw.org/install/windows.html
2. extract the zip. Copy the fftw3.h and libfftw3-3.dll and put them into your project folder, the folder same as *.pro file.
LIBS += "PWD/libfftw3-3.dll"
4. in QT creator, on the “projects” column, right click on the project name and select “Add Existing Files…”, choose the fftw.h
#include "fftw.h"
6. got to the build directory, usually named as “Build-…..”, you will see “debug” or “release” folders. copy the libfftw3-3.dll and paste it in the “debug” or “release” folder depends on you are debug or release.

The step 6 is very important, without so, the program will not even started, while the compilation is completely normal.

Electron spin resonance in power sample

The topic was discussed in a paper “Dynamic nuclear polarization by photoexcited-triplet electron spins in polycrystalline samples” by K. Takeda, K. Takegoshi, and T. Terao. I am going to give some supplementary on how the calculation can be done.

In a spin-1 system, the spin operator $S_x, S_y, S_z$ and the rotation operator are easy to calculate.

The zero-field splitting Hamiltonian in the ordinary basis of $S_z$ is

$\displaystyle H_{zfs} = \begin{pmatrix} D/3 & 0 & E \\ 0 & -2D/3 & 0 \\ E & 0 & D/3 \end{pmatrix}$

The eigen vectors in the paper are

$\displaystyle \left|X\right> = \frac{1}{\sqrt{2}}( \left|-\right> - \left|+\right> )$

$\displaystyle \left|Y\right> = \frac{i}{\sqrt{2}}( \left|-\right> + \left|+\right> )$

$\displaystyle \left|Z\right> = \left|0\right>$

with eigen values $X = D/3 -E$, $Y = D/3 +E$, and $Z = -2D/3$. The convention that the $\left| Y \right>$ has an imaginary number $i$ is not necessary. In the paper, there are typos on the eigen values.

In other way, the zero-field splitting Hamiltonian can be diagonalized into

$\displaystyle H_{zfs} = P\cdot H_{XYZ} \cdot P^{-1}$

$\displaystyle H_{XYZ} = \begin{pmatrix} X & 0 & 0 \\ 0 & Y & 0 \\ 0 & 0 & Z \end{pmatrix}$

$\displaystyle P = \frac{1}{\sqrt{2}} \begin{pmatrix} -1 & i & 0 \\ 0 & 0 & \sqrt{2} \\ 1 & i & 0 \end{pmatrix}$

where $P$ is the transformation changing from ordinary basis to eigen basis.

Given that the population of the eigen basis is $w = ( w_X, w_Y, w_Z)$, the density matrix is

$\displaystyle \rho_{XYZ} = \sum_{i = X, Y, Z} \left| i \right>\left< i \right| w_i$

The corresponding population in ordinary basis is the diagonal elements of

$\rho = P \cdot \rho_{XYZ} \cdot P^{-1}$

In the paper, it gives $w = (0.76, 0.16, 0.08)$. In order to get the population when the molecule X-axis is aligned with the Z-axis Lab frame, the $\rho$ has to be rotated by 90 degree, and we have the population $w' = ( 0.12, 0.76, 0.12)$.

The Zeeman Hamiltonian is $H_B = \omega S_Z$, where $\omega = \gamma_e B \hbar$. In the paper, the author fix the crystal and rotating the external magnetic field. To express the $H_B$ into eigen basis. We simply apply a transformation as the density matrix. But first we need to rotate the magnetic field. The rotation is simple, as the direction of the magnetic field is $(\theta, \phi)$. Thus, the Zeeman Hamiltonian is $H_B = \omega \hat{n} \cdot \vec{S}$.

The transformation is

$H_{B_XYZ} = P^{-1} \cdot H_B \cdot P$

again, the total Hamiltonian is difference from the author.

$\displaystyle H = \begin{pmatrix} X & -i \omega cos(\theta) & i \omega cos(\phi) sin(\theta) \\ i \omega cos(\theta) & Y & -i\omega sin(\phi) sin(\theta) \\ -i\omega cos(\phi) sin(\theta) & i \omega sin(\phi) sin(\theta) & Z \end{pmatrix}$

The difference is because the normal vector $\hat{n} = ( sin(\phi) sin(\theta), cos(\phi) sin(\theta), cos(\theta) )$ is used in the paper.

After that, the calculation for the power sample is the same.

Instead of using integration, we can generated the power using isotropic distribution. I generated 10k molecule, the result of 9.65 GHz microwave is

where the y axis is the magnetic field in Tesla. The result is similar to the paper. I think that their calculation somehow changed the eigen value, so that the eigen vectors were matching with the eigen energy. And the notation of the normal vector and the constant factor of the eigen vectors does not matter.

Phase shift of elastics scattering

I found that the derivation of most “google result” is not clear enough. So here is my derivation. Before process, people may need to review the pervious post.

Most people start from,

$u_l = A( \hat{h}_l^- - S_l \hat{h}_l^+ )$

that annoying me because it is somehow not “natural”. Why there is a “minus” sign? Why the $\hat{h}_l^-$ is the first term? For my self, a more natural way is,

$u_l = a \hat{h}_l^+ + b \hat{h}_l^-$

but that is not so natural at all, because in numerical calculation, for simplicity, there is no complex number, we only have,

$u_l = \alpha \hat{j}_l + \beta \hat{n}_l$

The first term is alway there as it is the free solution and bounded at $r = 0$. the second term is caused by the potential.

The goal is to find a solution take the form

$\displaystyle \psi = A \left( e^{i \vec{k} \cdot \vec{r}} + f(\theta) \frac{e^{ikr}}{r} \right)$

where the first term is free wave and the second term is scattered wave. The solution for elastics scattering is

$\displaystyle \psi = \sum C_l P_l (\cos\theta) \frac{u_l}{kr} \rightarrow \sum C_l P_l(\cos\theta) (\alpha \hat{j}_l + \beta \hat{n}_l)$

we used the substitution,

$\displaystyle R_l(r) = \frac{u_l(\rho)}{\rho}, \rho = kr$.

The radial function can be solved using Rungu-Kutta method on the equation,

$\displaystyle \frac{d^2}{d\rho^2} u_l = \frac{2 m_\mu}{\hbar^2} (V-E) u_l + \frac{l(l+1)}{\rho^2}$

and the solution of $u_l$ at far away is,

$u_l \rightarrow \alpha \hat{j}_l + \beta \hat{n}_l$.

the arrow means $r \rightarrow \infty$. So, the problem is how to rewrite the solution. In the way, we will see how the phase shift or the S-matrix was found.

The free solution is the spherical wave,

$\displaystyle e^{i \vec{k} \cdot \vec{r}} = \sum_l (2l+1) i^l P_l(\cos\theta) j_l(kr)$

The spherical Bessel function $j_l(kr)$ cna be express as Heankel function

$h_l^{\pm} = n_l \pm i j_l \rightarrow e^{\pm i (kr - l \frac{\pi}{2})}$

The $+$ sign is outgoing wave.

$\displaystyle u_l \rightarrow (\alpha \hat{j}_l + \beta \hat{n}_l)$

$\displaystyle = \frac{\alpha}{2i} (\hat{h}_l^{+} - \hat{h}_l^{-}) + \frac{\beta}{2}(\hat{h}_l^{+} + \hat{h}_l^{-})$

$\displaystyle = \frac{\alpha + i \beta}{2i} (\hat{h}_l^{+} - \hat{h}_l^{-}) + \beta \hat{h}_l^{+}$

$\displaystyle = (\alpha - i \beta ) \left( \frac{\hat{h}_l^+ - \hat{h}_l^-}{2i} + \frac{\beta}{\alpha - i \beta} \hat{h}_l^+\right)$

$\displaystyle = (\alpha - i \beta ) \left( \hat{j}_l + \frac{\beta}{\alpha - i \beta} \hat{h}_l^+\right)$

Since the $u_l$ should be normalized, we can se $\alpha = \cos \delta$ and $\beta = \sin\delta$.

$\displaystyle \frac{\beta}{\alpha - i \beta } = \sin(\delta) e^{i\delta}$

We put $u_l$ back

$\displaystyle \psi \rightarrow \sum_l C_l P_l (cos\theta)(\alpha - i \beta ) \left( \hat{j}_l + \sin(\delta) e^{i\delta} \hat{h}_l^+\right)$

By setting

$\displaystyle C_l = A i^l \frac{2l+1}{\alpha - i \beta}$,

we have the first term is the free wave function. In the second term, $\hat{h}_l^+ \rightarrow e^{i(kr - l \frac{\pi}{2}}) / kr$. Notice that

$e^{i l \frac{\pi}{2}} = i^{-l}$

That cancel the $i^l$ term in $C_l$. And we have

$\displaystyle f(\theta) = \sum (2l+1) P_l (\cos\theta) \frac{\sin(\delta) e^{i\delta}}{k}$

some people will write the $u_l$ as $\hat{h}_l^{\pm}$ and the S-matrix,

$\displaystyle u_l = \frac{\alpha + i \beta} {2i} \hat{h}_l^+ - \frac{\alpha - i \beta}{2i} \hat{h}_l^-$

$\displaystyle = -\frac{\alpha - i \beta}{2i} \left( \hat{h}_l^- - \frac{\alpha + i \beta}{\alpha - i \beta} \hat{h}_l^+ \right)$

$\displaystyle = A' (\hat{h}_l^- - S_l \hat{h}_l^+)$

where

$\displaystyle S_l =\frac{\alpha + i \beta}{\alpha - i \beta} = e^{2i\delta}$.

Remember that this is the S-matrix for elastics scattering.

Cross section IV

A more microscopic view of the cross section is

$\displaystyle \frac{dN}{dS} \sigma \rho dt dS= dn$,

where $N$ is number of incident particle, $\sigma$ is the total cross section, $\rho$ is the particle density in the target, $S$ is the area, $t$ is the thickness of the target, and $n$ is particle detected.

The cross section is a constant, so, after integration,

$\displaystyle \sigma = \frac{n}{N \rho t}$.

If the beam is not uniform, $dN/dS = f(S)$ is a function of $S$. The integration has to be careful.

A reminder of Exponential function

Suddenly found that exponential function is very wonderful. It is the eigen-function of differential operator. because of that, it appears in almost everywhere as differential operator constantly appears.

The definition of the exponential $e$ is,

$\displaystyle e = \lim_{n\to\infty} \left(1+\frac{1}{n} \right)^n$

this is the 1% growth continually in a unit of time. It is easy to see that

$\displaystyle e^x = \lim_{n\to\infty} \left(1+\frac{x}{n} \right)^n$

$e^x e^y = e^{x+y}$

expanding,

$\displaystyle \left(1+\frac{x}{n} \right)^n = \sum_{r=0}^{n} C_r^n \left(\frac{x}{n} \right)^r = \sum_{r=0}^{n} \left(\prod_{k=0}^{r-1} \frac{n-k}{n} \right)\frac{x^r}{r!}$,

Take $n\to\infty$,

$\displaystyle e^x = \sum_{r=0}^{\infty} \frac{x^r}{r!}$

It is easy to see that

$\displaystyle \frac{d}{dx} e^x = e^x$