11 The Hadamard Transform: Encoding and Decoding
Contents
11 The Hadamard Transform: Encoding and Decoding#
# import all python add-ons etc that will be needed later on
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from sympy import *
from scipy.integrate import quad
init_printing() # allows printing of SymPy results in typeset maths format
plt.rcParams.update({'font.size': 14}) # set font size for plots
Motivation and concept#
To obtain the average value of a quantity, such as weight, individual measurements are usually made, then added together and divided by the number of measurements made. However, this is not the only way of obtaining the average. You may be familiar with the method of weighing several objects at a time and perhaps know also that doing this will reduce the error in the average value. This multiple weighing method is an example of using a Hadamard transform and anything that can be measured in groups can be treated in the same way, for example, a spectrum or an image, thus the method is quite general. The reason that this multiple measuring method works is that the error is introduced by the balance not by the objects being weighed. A large weight therefore has the same error associated with its measurement as a smaller one does.
The reason for doing any transform experiment is always the same and either this is to achieve an improvement in signal to noise, or, a reduction in the time taken to do an experiment at a fixed signal to noise, which is effectively the same thing. Normally for
11.1 The Hadamard Transform#
The Hadamard transform is a purely discrete transform and instead of forward and back transforming, as in the Fourier transform, the equivalent steps are encoding and decoding. The encoding is done by adding several measurements together according to a set of rules or algorithm. The rule is always written down as a matrix, two forms of which
11.2 Encoding and decoding with matrices#
Suppose that there are three samples to be weighed of masses
These equations can be solved simultaneously to find the weights (
Note how the pattern is the same as that of the coefficients in the equations and that it cycles around so that each column is related to the next by one position of cyclic rotation. As a matrix equation, equation (47) is
and to solve for (column) vector
To show that this works, suppose that the
Of course if there are more than three values then a different
# check on matrix inverse
S = symbols('S')
S = Matrix([[1, 1, 0], [1, 0, 1], [0, 1,1]])
S**(-1)
11.3 Signal to Noise improvement#
To see why this method works to improve signal to noise some error has to be added to each measurement. If each measurement has a standard deviation determined by the instrument (scales) used then a measurement of
which is used to estimate the overall error. If
z,m1,m2,m3,sigma = symbols('z m1 m2 m3 sigma')
z = Matrix([m1 + m2 + sigma, m1 + m3 + sigma,m2 + m3 + sigma]) # define matrix
psi = S**(-1)*z
psi
subtracting each
11.4 Implementation#
Instead of weights, suppose that a spectrum is to be measured. To do the experiment, the detector is placed at the focusing plane of a spectrometer, the exit slits of which are removed and then a mask consisting of strips of opaque (0) and transparent (1) regions is placed there instead, see Fig. 49. At the first position the total amount of light falling on the detector is measured, this is
Figure 49. The pattern of the mask replaces the slits of the spectrometer. The detector measures all the light transmitted by the mask at each position. Each mask is rotated by one element from the previous one. All the mask elements when placed together form the
In Fig. 50 is shown a simulated comparison of data taken in the normal way and with that taken using the Hadamard encoding method. The noise on the detector is normally distributed with a mean of
The Hadamard technique has recently been applied by the author and colleagues to time-resolved x-ray diffraction (Nat. Methods. 2014 Nov; 11(11): 1131-1134. ) and has also been applied to time-resolved spectroscopic measurements (Appl. Spectros. 2016,70,1292-1299). To perform a time-resolved absorption measurement, a train of pulses whose intensity is in the
The timescale of the measurement is set by the spacing between any two pulses in the sequence and the total measurement time is
Figure 50. Comparison of Hadamard and normal (conventional) on the same set of data. The initial ideal data is shown in both figures as the blue line, the ‘measured’ data is in red.
11.5 Constructing the matrix#
Harwit & Sloane (1979) give several methods by which to construct the
The Quadratic residue
A flow diagram to make a row of the
# check for valid sequence length. First few values are 3, 7, 11, 19, 23, 31, 43, 47,
#------------
def valid_seq_length(n):
maxi = 200
Hseq = np.zeros(maxi,dtype=int)
for i in range(maxi): # produce Hadamard sequence numbers
for m in range(0,maxi):
if isprime(i) and i == 4*m + 3:
Hseq[i]=i
pass
if n in Hseq[0:maxi]:
is_ok = True
else:
is_ok = False
return is_ok
#------------
def isprime(n): # check if integer n is a prime, range starts with 2 and only needs to go up the squareroot of n
for x in range(2, int(n**0.5)+1):
if n % x == 0:
return False
return True
#------------
def make_S_mat(Srow):
n = len(Srow)
S = np.zeros((n,n),dtype = int)
for i in range(n):
for j in range(n):
S[i,j] = Srow[n-1-j]
Srow = np.roll(Srow, -1) # rotate by 1 element at a time
pass
return S
#------------
def quadratic_hadamard(n): # quadratic residue method, this generates an S matrix.
alist = np.zeros(n,dtype=int)
for i in range(0,(n-1)//2): # integer division
alist[(i+1)*(i+1) % n] = 1.0 # alist = hadamard Srow need only go to half range of n as indices are symmetric
alist[0] = 1.0
Srow = list(alist)
S = make_S_mat(Srow[::-1]) # ***** reverse to get 1's at start; no other reason ****
return S # returns S matrix
#------------
print('Hadamard S matrices by Quadratic residue method')
for i in range(1,12): # print out S matrices up to size of 12 x 12
if valid_seq_length(i):
S = quadratic_hadamard(i)
print(i)
if i <= 32 : # for larger matrices print just first line
print('\n'.join( [''.join(['{:2}'.format(item) for item in row] ) for row in S] ) )
else:
print(''.join(['{:2}'.format(item) for item in S[0]]))
xs=''.join( str(S[0][i]) for i in range(len(S[0])) )
print(hex(int(xs,2)))
pass
Hadamard S matrices by Quadratic residue method
3
1 1 0
0 1 1
1 0 1
7
1 1 1 0 1 0 0
0 1 1 1 0 1 0
0 0 1 1 1 0 1
1 0 0 1 1 1 0
0 1 0 0 1 1 1
1 0 1 0 0 1 1
1 1 0 1 0 0 1
11
1 1 0 1 1 1 0 0 0 1 0
0 1 1 0 1 1 1 0 0 0 1
1 0 1 1 0 1 1 1 0 0 0
0 1 0 1 1 0 1 1 1 0 0
0 0 1 0 1 1 0 1 1 1 0
0 0 0 1 0 1 1 0 1 1 1
1 0 0 0 1 0 1 1 0 1 1
1 1 0 0 0 1 0 1 1 0 1
1 1 1 0 0 0 1 0 1 1 0
0 1 1 1 0 0 0 1 0 1 1
1 0 1 1 1 0 0 0 1 0 1
The figure below shows two
In the calculation the inverse of the
where