CSC 726, Parallel Algorithms
QR Factorization
This page contains a list of examples and helpful code
fragments for running your QR factorization code.
Recall that our linear system is:
g = H f
where f is a vector size n*n by 1,
H is an matrix size 2*n*n by n*n, and
g is a vector size 2*n*n by 1. Note that
the system is over determined, so you will need
to compute leastsquares solution. QR factoriziation
is an important method for least squares problems
such as this one. A summary of the method and
how it applies to least squares problems is
linked here.
The vector f contains the pixel values of an n by n,
image stacked by columns. Similarly, g contains the
stacked pixel values for two images. The vector g
represents two observations, and the matrix H represents
the blurring operator for both observations. The inverse problem
posed here is to find f. For the following files, use
your code to compute an estimate of f using g
and H. In the example files below, n = 16.
For purposes of code validation, the original
vector f is also provided. Your solution
should be very close to the given f. Later, you
will have the opportunity to solve problems for which
f is not available.
The following ".dat" files contain littleendian binary data.
Note on file format: These are binary files. The first 4 bytes
form a littleendian integer M, representing the number of rows in
the object. The second 4 bytes form a littleendian integer N,
representing the number of columns in the object.
The next M * N * sizeof(double) bytes are littleendian double
precision values in columnmajor order. A sample seqential
C function to read this format is
included here. Sample sequential C++
code is here. Sample C code to
write this format is included here.
Sample C++ code to write this format is here.
Also, you will want to convert matricies to a graphics format
for viewing. A sequential program named dat2ppm to
convert to the PPM graphics format is
compiled and (currently) installed in /var/tmp/bin
on bally.cs.wfu.edu. The program dat2ppm takes
two options for controlling the colormap: dat2ppm gray filename.dat
produces a grayscale image; dat2ppm jet filename.dat produces
a pseudocolor image using a colormap similar to Matlab's "jet" colormap.
You can view PPM
graphics files on bally with the display command
from the ImageMagick software package (already installed on bally).
Sample Input Data for QR Factorization
A small problem for getting things tested and working.
 g_le.dat , a simulated observation (512 x 1)
 h_le.dat , a simulated blurring operator (512 x 256)
 f_le.dat , the unknown. (256 x 1)
What do these images look like ? The vector f is
unstacked into a square 16 x 16 image. The vector g
is split into the top and bottom regions g1 and g2,
and each is unstacked into a square 16 x 16 image.
For illustration purposes, each image is enlarged
(via block replication). Then, we have:
The matrix h is illustrated in the 'jet'
colormap below. Blues and cool tones indicate low
values; oranges, reds and warm tones indicate high
values. The top and bottom halves of h are
structured: they are block circulant with
circulant blocks.

Matrix h 
A somewhat larger problem.
Now that your code is working smoothly, try a slightly larger problem, based on an image size 32 x 32.
For the 2048 x 1024 matrix "h_32.dat" above, what does the
QR factorization look like ?
A more interesting problem.
You can also try a larger problem, based on an image size 64 x 64.
An illconditioned, spatially varying problem.
Finally, we have arrived at our first real problem.
The following files present an illconditioned simulated
spatiallyvarying blur on 64x64 images. You will need to
include regularization to produce an acceptable restoration.
Sample C code to read the data file and to assist with Tikhonov
regularization is found in tikhonov.h
and in tikhonov.c. These functions
should also compile easily if you are using the C++ compiler,
with a minor modification to the syntax for the #include files.
Images: