Docstoc

IMAGE FUSION

Document Sample
IMAGE FUSION Powered By Docstoc
					                              ABSTRACT

      The objective of the image fusion is to combine the source images of

the same scene to form one composite image contains a more accurate

description. Image fusion methods can be classified as ― spatial domain

fusion ‖ and ― transform domain fusion ”. Averaging,        Brovey method,

principal component analysis (PCA) and IHS methods fall under spatial

domain approaches. We propose a multiresolution data fusion scheme based

on the PCA transform and the pixel-level weights wavelet transform. A linear

local mapping based on the PCA is used to create a new "origin" image of the

image fusion. Daubechies wavelet is chosen as wavelet basis. Wavelet based

fusion techniques have been reasonably effective in combining perceptually

important image features. Shift invariance wavelet transform is important in

ensuring robust sub band fusion. Applications of image fusion are Satellite

imaging, Medical imaging, Robot vision, Concealed weapon detection, Multi-

focus image fusion, Digital camera application, Battle field      monitoring.

      The aim of image fusion is to create new images that are more suitable

for the purposes of human visual perception, object detection and target

recognition.

      Experimental results confirm that the proposed algorithm is the best

image sharpening method and can best maintain the spectral information of

the original image. Also, the proposed technique performs better than the

other ones, more robust and effective, from both subjective visual effects

and objective statistical analysis results. The performance of the image

fusion is evaluated by normalized least square error, entropy, overall cross

entropy, standard deviation and mutual information.
 1 INTRODUCTION TO DIGITAL IMAGE PROCESSING

1.1 IMAGE:

              A digital image is a computer file that contains graphical

information instead of text or a program. Pixels are the basic building blocks

of all digital images. Pixels are small adjoining squares in a matrix across the

length and width of your digital image. They are so small that you don‘t see

the actual pixels when the image is on your computer monitor.

        Pixels are monochromatic. Each pixel is a single solid color that is

blended from some combination of the 3 primary colors of Red, Green, and

Blue. So, every pixel has a RED component, a GREEN component and BLUE

component. The physical dimensions of a digital image are measured in

pixels and commonly called pixel or image resolution. Pixels are scalable to

different physical sizes on your computer monitor or on a photo print.

However, all of the pixels in any particular digital image are the same size.

Pixels as represented in a printed photo become round slightly overlapping

dots.
      1.1.1 Pixel Values: As shown in this bitonal image, each pixel is

assigned a tonal value, in this example 0 for black and 1 for white.


      1.1.2 Pixel Dimensions: They are the horizontal and vertical

measurements of an image expressed in pixels. The pixel dimensions may be

determined by multiplying both the width and the height by the dpi. A digital

camera will also have pixel dimensions, expressed as the number of pixels

horizontally and vertically that define its resolution (e.g., 2,048 by 3,072).

Calculate the dpi achieved by dividing a document's dimension into the

corresponding pixel dimension against which it is aligned.




      Example:




Fig: An 8" x 10" document that is scanned at 300 dpi has the pixel

dimensions of 2,400 pixels (8" x 300 dpi) by 3,000 pixels (10" x 300 dpi).
1.2 IMAGES IN MATLAB:




       The basic data structure in MATLAB is the array, an ordered set of real

or complex elements. This object is naturally suited to the representation of

images, real-valued ordered sets of color or intensity data.

       MATLAB stores most images as two-dimensional arrays (i.e., matrices),

in which each element of the matrix corresponds to a single pixel in the

displayed image. (Pixel is derived from picture element and usually denotes a

single dot on a computer display.)

       For example, an image composed of 200 rows and 300 columns of

different colored dots would be stored in MATLAB as a 200-by-300 matrix.

Some images, such as color images, require a three-dimensional array,

where the first plane in the third dimension represents the red pixel intensities,

the second plane represents the green pixel intensities, and the third plane

represents the blue pixel intensities. This convention makes working with

images in MATLAB similar to working with any other type of matrix data, and

makes the full power of MATLAB available for image processing applications.
1.3 IMAGE REPRESENTATION:


   An image is stored as a matrix using standard Matlab matrix conventions.

There are four basic types of images supported by Matlab:


   1.3.1 Binary images

   1.3.2 Intensity images

   1.3.3 RGB images

   1.3.4 Indexed images


1.3.1 Binary Images:


In a binary image, each pixel assumes one of only two discrete values: 1 or 0.

A binary image is stored as a logical array. By convention, this documentation

uses the variable name BW to refer to binary images.


The following figure shows a binary image with a close-up view of some of the

pixel values.




                    Fig: Pixel Values in a Binary Image
1.3.2 Intensity images:


      A grayscale image (also called gray-scale, gray scale, or gray-level) is

a data matrix whose values represent intensities within some range. MATLAB

stores a grayscale image as an individual matrix, with each element of the

matrix corresponding to one image pixel. By convention, this documentation

uses the variable name I to refer to grayscale images.


      The matrix can be of class uint8, uint16, int16, single, or double. While

grayscale images are rarely saved with a color map, MATLAB uses a color

map to display them.


      For a matrix of class single or double, using the default grayscale color

map, the intensity 0 represents black and the intensity 1 represents white. For

a matrix of type uint8, uint16, or int16, the intensity intmin (class (I))

represents black and the intensity intmax (class (I)) represents white.


The figure below depicts a grayscale image of class double.




   Fig: Pixel Values in a Grayscale Image Define Gray Levels
1.3.3 RGB Images:


       A color image is an image in which each pixel is specified by three

values — one each for the red, blue, and green components of the pixel's

color. MATLAB store color images as an m-by-n-by-3 data array that defines

red, green, and blue color components for each individual pixel. Color images

do not use a color map. The color of each pixel is determined by the

combination of the red, green, and blue intensities stored in each color plane

at the pixel's location.


       Graphics file formats store color images as 24-bit images, where the

red, green, and blue components are 8 bits each. This yields a potential of 16

million colors. The precision with which a real-life image can be replicated has

led to the commonly used term color image.


       A color array can be of class uint8, uint16, single, or double. In a color

array of class single or double, each color component is a value between 0

and 1. A pixel whose color components are (0, 0, 0) is displayed as black, and

a pixel whose color components are (1, 1, 1) is displayed as white. The three

color components for each pixel are stored along the third dimension of the

data array. For example, the red, green, and blue color components of the

pixel (10,5) are stored in RGB(10,5,1), RGB(10,5,2), and RGB(10,5,3),

respectively.


The following figure depicts a color image of class double.
Fig: Color Planes of a True color Image
1.3.4 Indexed Images:


       An indexed image consists of an array and a color map matrix. The

pixel values in the array are direct indices into a color map. By convention,

this documentation uses the variable name X to refer to the array and map to

refer to the color map.


       The color map matrix is an m-by-3 array of class double containing

floating-point values in the range [0, 1]. Each row of map specifies the red,

green, and blue components of a single color. An indexed image uses direct

mapping of pixel values to color map values. The color of each image pixel is

determined by using the corresponding value of X as an index into map.


       A color map is often stored with an indexed image and is automatically

loaded with the image when you use the imread function. After you read the

image and the color map into the MATLAB workspace as separate variables,

you must keep track of the association between the image and color map.

However, you are not limited to using the default color map--you can use any

color map that you choose.


       The relationship between the values in the image matrix and the color

map depends on the class of the image matrix. If the image matrix is of class

single or double, it normally contains integer values 1 through p, where p is

the length of the color map. The value 1 points to the first row in the color map,

the value 2 points to the second row, and so on. If the image matrix is of class
       Logical, uint8 or uint16, the value 0 points to the first row in the color

map, the value 1 points to the second row, and so on.


The following figure illustrates the structure of an indexed image. In the figure,

the image matrix is of class double, so the value 5 points to the fifth row of the

color map.




    Fig: Pixel Values Index to Color map Entries in Indexed

                                    Images
1.4 DIGITAL IMAGE FILE TYPES:




The 5 most common digital image file types are as follows:




1.4.1 JPEG: It is a compressed file format that supports 24 bit color

(millions of colors). This is the best format for photographs to be shown on the

web or as email attachments. This is because the color informational bits in

the computer file are compressed (reduced) and download times are

minimized.




1.4.2 GIF: It is an uncompressed file format that supports only 256 distinct

colors. Best used with web clip art and logo type images. GIF is not suitable

for photographs because of its limited color support.




1.4.3 TIFF: It is an uncompressed file format with 24 or 48 bit color support.

Uncompressed means that all of the color information from your scanner or

digital camera for each individual pixel is preserved when you save as TIFF.

TIFF is the best format for saving digital images that you will want to print. Tiff

supports embedded file information, including exact color space, output profile
information and EXIF data. There is a lossless compression for TIFF called

LZW. LZW is much like 'zipping' the image file because there is no quality loss.

An LZW TIFF decompresses (opens) with all of the original pixel information

unaltered.




1.4.4 BMP: It is a Windows (only) operating system uncompressed file

format that supports 24 bit color. BMP does not support embedded

information like EXIF, calibrated color space and output profiles. Avoid using

BMP for photographs because it produces approximately the same file sizes

as TIFF without any of the advantages of TIFF.




1.4.5 Camera RAW: is a lossless compressed file format that is

proprietary for each digital camera manufacturer and model. A camera RAW

file contains the 'raw' data from the camera's imaging sensor. Some image

editing programs have their own version of RAW too. However, camera RAW

is the most common type of RAW file. The advantage of camera RAW is that

it contains the full range of color information from the sensor. This means the

RAW file contains 12 to 14 bits of color information for each pixel. If you shoot

JPEG, you only get 8 bits of color for each pixel. These extra color bits make

shooting camera RAW much like shooting negative film. You have a little

more latitude in setting your exposure and a slightly wider dynamic range.
1.5 IMAGE COORDINATE SYSTEMS:




1.5.1 Pixel Coordinates:


       Generally, the most convenient method for expressing locations in an

image is to use pixel coordinates. In this coordinate system, the image is

treated as a grid of discrete elements, ordered from top to bottom and left to

right, as illustrated by the following figure.




                    Fig: The Pixel Coordinate System


       For pixel coordinates, the first component r (the row) increases

downward, while the second component c (the column) increases to the right.

Pixel coordinates are integer values and range between 1 and the length of

the row or column.


     There is a one-to-one correspondence between pixel coordinates and the

coordinates MATLAB uses for matrix subscripting. This correspondence

makes the relationship between an image's data matrix and the way the
image is displayed easy to understand. For example, the data for the pixel in

the fifth row, second column is stored in the matrix element (5, 2). You use

normal MATLAB matrix subscripting to access values of individual pixels.




     For example, the MATLAB code

    I (2, 15)


Returns the value of the pixel at row 2, column 15 of the image I.


1.5.2 Spatial Coordinates:


       In the pixel coordinate system, a pixel is treated as a discrete unit,

uniquely identified by a single coordinate pair, such as (5, 2). From this

perspective, a location such as (5.3, 2.2) is not meaningful.


       At times, however, it is useful to think of a pixel as a square patch.

From this perspective, a location such as (5.3, 2.2) is meaningful, and is

distinct from (5, 2). In this spatial coordinate system, locations in an image are

positions on a plane, and they are described in terms of x and y (not r and c

as in the pixel coordinate system).


The following figure illustrates the spatial coordinate system used for images.

Notice that y increases downward.
                  Fig: The Spatial Coordinate System




       This spatial coordinate system corresponds closely to the pixel

coordinate system in many ways. For example, the spatial coordinates of the

center point of any pixel are identical to the pixel coordinates for that pixel.


       There are some important differences, however. In pixel coordinates,

the upper left corner of an image is (1,1), while in spatial coordinates, this

location by default is (0.5,0.5). This difference is due to the pixel coordinate

system's being discrete, while the spatial coordinate system is continuous.

Also, the upper left corner is always (1,1) in pixel coordinates, but you can

specify a non default origin for the spatial coordinate system.


       Another potentially confusing difference is largely a matter of

convention: the order of the horizontal and vertical components is reversed in

the notation for these two systems. As mentioned earlier, pixel coordinates

are expressed as (r, c), while spatial coordinates are expressed as (x, y). In

the reference pages, when the syntax for a function uses r and c, it refers to

the pixel coordinate system. When the syntax uses x and y, it refers to the

spatial coordinate system.
1.6 DIGITAL IMAGE PROCESSING:


       Digital image processing is the use of computer algorithms to perform

image processing on digital images. As a subfield of digital signal processing,

digital image processing has many advantages over analog image

processing; it allows a much wider range of algorithms to be applied to the

input data, and can avoid problems such as the build-up of noise and signal

distortion during processing.


1.6.1 Image digitization:

       An image captured by a sensor is expressed as a continuous function

f(x,y) of two co-ordinates in the plane. Image digitization means that the

function f(x,y) is sampled into a matrix with M rows and N columns. The image

quantization assigns to each continuous sample an integer value. The

continuous range of the image function f(x,y) is split into K intervals. The finer

the sampling (i.e., the larger M and N) and quantization (the larger K) the

better the approximation of the continuous image function f(x,y).


1.6.2 Image Pre-processing:

   Pre-processing is a common name for operations with images at the

lowest level of abstraction -- both input and output are intensity images. These

iconic images are of the same kind as the original data captured by the sensor,

with an intensity image usually represented by a matrix of image function

values (brightness). The aim of pre-processing is an improvement of the

image data that suppresses unwanted distortions or enhances some image

features important for further processing. Four categories of image pre-
processing methods according to the size of the pixel neighborhood that is

used for the calculation of new pixel brightness:


          o   Pixel brightness transformations.

          o   Geometric transformations.

          o   Pre-processing methods that use a local neighborhood of the

              processed pixel.

          o   Image restoration that requires knowledge about the entire

              image.


1.6.3 Image Segmentation:


   Image segmentation is one of the most important steps leading to the

analysis of processed image data. Its main goal is to divide an image into

parts that have a strong correlation with objects or areas of the real world

contained in the image. Two kinds of segmentation


   1.    Complete segmentation: This results in set of disjoint regions

         uniquely corresponding with objects in the input image. Cooperation

         with higher processing levels which use specific knowledge of the

         problem domain is necessary.

   2.    Partial segmentation: in which regions do not correspond directly

         with image objects. Image is divided into separate regions that are

         homogeneous with respect to a chosen property such as brightness,

         color, reflectivity, texture, etc. In a complex scene, a set of possibly

         overlapping   homogeneous      regions may result.       The   partially

         segmented image must then be subjected to further processing, and
         the final image segmentation may be found with the help of higher

         level information.


Segmentation methods can be divided into three groups according to the

dominant features they employ


          1. First is global knowledge about an image or its part; the

             knowledge is usually represented by a histogram of image

             features.

          2. Edge-based segmentations form the second group; and

          3. Region-based segmentations


1.6.4 Image enhancement:


      The aim of image enhancement is to improve the interpretability or

perception of information in images for human viewers, or to provide `better'

input for other automated image processing techniques. Image enhancement

techniques can be divided into two broad categories:


1. Spatial domain methods, which operate directly on pixels, and

2. Frequency domain methods, which operate on the Fourier transform of an

image.


Unfortunately, there is no general theory for determining what `good‘ image

enhancement is when it comes to human perception. If it looks good, it is

good! However, when image enhancement techniques are used as pre-

processing tools for other image processing techniques, then quantitative

measures can determine which techniques are most appropriate.
            2 INTRODUCTION TO IMAGE FUSION:



2.1 INTRODUCTION:



       In computer vision, Multisensor Image Fusion is the process of

combining relevant information from two or more images into a single image.

The resulting image will be more informative than any of the input images. In

remote sensing applications, the increasing availability of space borne

sensors gives a motivation for different image fusion algorithms. Several

situations in image processing require high spatial and high spectral

resolution in a single image. Most of the available equipment is not capable of

providing such data convincingly. The image fusion techniques allow the

integration of different information sources. The fused image can have

complementary spatial and spectral resolution characteristics. But, the

standard image fusion techniques can distort the spectral information of the

multispectral data, while merging.


       In satellite imaging, two types of images are available. The

panchromatic image acquired by satellites is transmitted with the maximum

resolution available and the multispectral data are transmitted with coarser

resolution. This will be usually, two or four times lower. At the receiver station,

the panchromatic image is merged with the multispectral data to convey more

information.
Many methods exist to perform image fusion. The very basic one is the high

pass filtering technique. Later techniques are based on DWT, uniform rational

filter bank, and Laplacian pyramid.


       Multisensor data fusion has become a discipline to which more and

more general formal solutions to a number of application cases are demanded.

Several situations in image processing simultaneously require high spatial

and high spectral information in a single image. This is important in remote

sensing. However, the instruments are not capable of providing such

information either by design or because of observational constraints. One

possible solution for this is data fusion

2.2 STANDARD IMAGE FUSION METHODS:


   Image fusion methods can be broadly classified into two - spatial domain

fusion and transform domain fusion. The fusion methods such as averaging,

Brovey method, principal component analysis (PCA) and IHS based methods

fall under spatial domain approaches. Another important spatial domain fusion

method is the high pass filtering based technique. Here the high frequency

details are injected into up sampled version of MS images. The disadvantage

of spatial domain approaches is that they produce spatial distortion in the

fused image. Spectral distortion becomes a negative factor while we go for

further processing, such as classification problem, of the fused image. The

spatial distortion can be very well handled by transform domain approaches

on image fusion.


   The multiresolution analysis has become a very useful tool for analyzing

remote sensing images. The discrete wavelet transform has become a very
useful tool for fusion. Some other fusion methods are also there, such as

Laplacian pyramid based, Curve let transform based etc. These methods

show a better performance in spatial and spectral quality of the fused image

compared        to      other        spatial     methods        of      fusion.

2.2.1 Applications:


   1. Image Classification

   2. Aerial and Satellite imaging

   3. Medical imaging

   4. Robot vision

   5. Concealed weapon detection

   6. Multi-focus image fusion

   7. Digital camera application

   8. Concealed weapon detection

   9. Battle field monitoring


2.3 SATELLITE IMAGE FUSION:


Several methods are there for merging satellite satilite images. In satellite

imagery we can have two types of images


      Panchromatic images - An image collected in the broad visual

       wavelength range but rendered in black and white.

      Multispectral images - Images optically acquired in more than one

       spectral or wavelength interval. Each individual image is usually of the

       same physical area and scale but of a different spectral band.
The SPOT PAN satellite provides high resolution (10m pixel) panchromatic

data while the LANDSAT TM satellite provides low resolution (30m pixel)

multispectral images. Image fusion attempts to merge these images and

produce a single high resolution multispectral image.




The standard merging methods of image fusion are based on Red-Green-

Blue (RGB) to Intensity-Hue-Saturation (IHS) transformation. The usual steps

involved in satellite image fusion are as follows:




   1. Register the low resolution multispectral images to the same size as

       the panchromatic image

   2. Transform the R,G and B bands of the multispectral image into IHS

       components

   3. Modify the panchromatic image with respect to the multispectral image.

       This is usually performed by Histogram Matching of the panchromatic

       image with Intensity component of the multispectral images as

       reference

   4. Replace the intensity component by the panchromatic image and

       perform inverse transformation to obtain a high resolution multispectral

       image.
2.4 MEDICAL IMAGE FUSION:




      Image fusion has recently become a common term used within medical

diagnostics and treatment. The term is used when patient images in different

data formats are fused. These forms can include magnetic resonance image

(MRI), computed tomography (CT), and positron emission tomography (PET).

In radiology and radiation oncology, these images serve different purposes.

For example, CT images are used more often to ascertain differences in

tissue density while MRI images are typically used to diagnose brain tumors.




      For accurate diagnoses, radiologists must integrate information from

multiple image formats. Fused, anatomically-consistent images are especially

beneficial in diagnosing and treating cancer. Companies such as Keosys,

MIMvista, IKOE, and Brain LAB have recently created image fusion software

to use in conjunction with radiation treatment planning systems. With the

advent of these new technologies, radiation oncologists can take full

advantage of intensity modulated radiation therapy (IMRT). Being able to

overlay diagnostic images onto radiation planning images results in more

accurate IMRT target tumor volumes.
2.5 IMAGE FUSION ALGORITHMS:



      The details of wavelets and PCA algorithm and their use in image

fusion along with simple average fusion algorithm are described in this section.



2.5.1 Principal Component Analysis:

      The PCA involves a mathematical procedure that transforms a number

of correlated variables into a number of uncorrelated variables called principal

components. It computes a compact and optimal description of the data set.

The first principal component accounts for as much of the variance in the data

as possible and each succeeding component accounts for as much of the

remaining variance as possible. First principal component is taken to be along

the direction with the maximum variance. The second principal component is

constrained to lie in the subspace perpendicular of the first. Within this

subspace, this component points the direction of maximum variance. The third

principal component is taken in the maximum variance direction in the

subspace perpendicular to the first two and so on. The PCA is also called as

Karhunen-Loève transform or the Hotelling transform. The PCA does not have

a fixed set of basis vectors like FFT, DCT and wavelet etc. and its basis

vectors depend on the data set.

Let X be a d-dimensional random vector and assume it to have zero empirical

mean. Orthonormal projection matrix V would be such that Y=V TX with the
following constraints. The covariance of Y, i.e., cov(Y) is a diagonal and

inverse of V is equivalent to its transpose ( V-1 = VT).




Using matrix algebra




                                                   …………… (1)

Multiplying both sides of above eqn by V, one gets


One could write V as V= [V1, V2… , Vd] and

cov(Y) as




                                                  …………… (2)

Substituting Eqn (1) into the Eqn (2) gives

[λ1V1, λ2V2, …. , λdVd] = [cov(X) V1, cov(X) V2, … , cov(X) Vd] ………. (3)

This could be rewritten as

                           λiVi = cov(X) V1 ……… (4)

where i=1,2,...,d and Vi is an eigenvector of cov(X ) .
2.5.2 PCA Algorithm:

Let the source images (images to be fused) be arranged in two-column

vectors. The steps followed to project this data into 2-D subspaces are:

1. Organize the data into column vectors. The resulting matrix Z is of

dimension 2 x n.

2. Compute the empirical mean along each column. The empirical mean

vector M has a dimension of 1 x 2.

3. Subtract the empirical mean vector M from each column of the data matrix

Z. The resulting matrix X is of dimension 2 x n.

4. Find the covariance matrix C of X i.e. =XXT mean of expectation = cov(X)

5. Compute the eigenvectors V and eigenvalue D of C and sort them by

decreasing eigenvalue. Both V and D are of dimension 2 x 2.

6. Consider the first column of V which corresponds to larger eigenvalue to

compute P1 and P2 as




And
2.5.3 Image Fusion by PCA:

      The information flow diagram of PCA-based image fusion algorithm is

shown in figure below. The input images (images to be fused) I1 (x, y) and I2

(x, y) are arranged in two column vectors and their empirical means are

subtracted. The resulting vector has a dimension of n x 2, where n is length of

the each image vector. Compute the eigenvector and eigenvalues for this

resulting vector are computed and the eigenvectors corresponding to the

larger eigenvalue obtained. The normalized components P1 and P2 (i.e., P1 +

P2 = 1) using equation (3) are computed from the obtained eigenvector. The

fused image is:




   Figure: Information flow diagram in image fusion scheme

                             employing PCA.

Image Fusion by Simple Average:

This technique is a basic and straightforward technique and fusion could be

achieved by simple averaging corresponding pixels in each input image as:
2.6 IMAGE FUSION BY WAVELET TRANSFORMS:


2.6.1 Fourier analysis:


      Signal analysts already have at their disposal an impressive arsenal of
tools. Perhaps the most well-known of these is Fourier analysis, which breaks
down a signal into constituent sinusoids of different frequencies. Another way
to think of Fourier analysis is as a mathematical technique for transforming
our view of the signal from time-based to frequency-based.




                                     Figure 2
     For many signals, Fourier analysis is extremely useful because the
signal‘s frequency content is of great importance. So why do we need other
techniques, like wavelet analysis?


     Fourier analysis has a serious drawback. In transforming to the
frequency domain, time information is lost. When looking at a Fourier
transform of a signal, it is impossible to tell when a particular event took place.
If the signal properties do not change much over time — that is, if it is what is
called a stationary signal—this drawback isn‘t very important. However, most
interesting   signals   contain   numerous      non     stationary   or   transitory
characteristics: drift, trends, abrupt changes, and beginnings and ends of
events. These characteristics are often the most important part of the signal,
and Fourier analysis is not suited to detecting them.
2.6.2 Short-Time Fourier Analysis:


     In an effort to correct this deficiency, Dennis Gabor (1946) adapted the
Fourier transform to analyze only a small section of the signal at a time—a
technique called windowing the signal. Gabor‘s adaptation, called the Short-
Time Fourier Transform (STFT), maps a signal into a two-dimensional
function of time and
frequency.




                                  Figure 3


      The STFT represents a sort of compromise between the time- and
frequency-based views of a signal. It provides some information about both
when and at what frequencies a signal event occurs. However, you can only
obtain this information with limited precision, and that precision is determined
by the size of the window. While the STFT compromise between time and
frequency information can be useful, the drawback is that once you choose a
particular size for the time window, that window is the same for all frequencies.
Many signals require a more flexible approach—one where we can vary the
window size to determine more accurately either time or frequency.
2.7 WAVELET ANALYSIS:


     Wavelet analysis represents the next logical step: a windowing
technique with variable-sized regions. Wavelet analysis allows the use of long
time intervals where we want more precise low-frequency information, and
shorter regions where we want high-frequency information.




                                  Figure 4
Here‘s what this looks like in contrast with the time-based, frequency-based,
and STFT views of a signal:




                                  Figure 5
You may have noticed that wavelet analysis does not use a time-frequency
region, but rather a time-scale region. For more information about the concept
of scale and the link between scale and frequency, see ―How to Connect
Scale to Frequency?‖
2.7.1 What Can Wavelet Analysis Do?


     One major advantage afforded by wavelets is the ability to perform local
analysis, that is, to analyze a localized area of a larger signal. Consider a
sinusoidal signal with a small discontinuity — one so tiny as to be barely
visible. Such a signal easily could be generated in the real world, perhaps by
a power fluctuation or a noisy switch.




                                   Figure 6
   A plot of the Fourier coefficients (as provided by the fft command) of this
signal shows nothing particularly interesting: a flat spectrum with two peaks
representing a single frequency. However, a plot of wavelet coefficients
clearly shows the exact location in time of the discontinuity.




                                   Figure 7
                          Wavelet analysis is capable of revealing aspects of
data that other signal analysis techniques miss, aspects like trends,
breakdown points, discontinuities in higher derivatives, and self-similarity.
Furthermore, because it affords a different view of data than those presented
by traditional techniques, wavelet analysis can often compress or de-noise a
signal without appreciable degradation. Indeed, in their brief history within the
signal processing field, wavelets have already proven themselves to be an
indispensable addition to the analyst‘s collection of tools and continue to enjoy
a burgeoning popularity today.




2.7.2 What Is Wavelet Analysis?



     Now that we know some situations when wavelet analysis is useful, it is
worthwhile asking ―What is wavelet analysis?‖ and even more fundamentally,
―What is a wavelet?‖
A wavelet is a waveform of effectively limited duration that has an average
value of zero.
Compare wavelets with sine waves, which are the basis of Fourier analysis.
Sinusoids do not have limited duration — they extend from minus to plus
infinity. And where sinusoids are smooth and predictable, wavelets tend to be
irregular and asymmetric.




                                   Figure 8
     Fourier analysis consists of breaking up a signal into sine waves of
various frequencies. Similarly, wavelet analysis is the breaking up of a signal
into shifted and scaled versions of the original (or mother) wavelet. Just
looking at pictures of wavelets and sine waves, you can see intuitively that
signals with sharp changes might be better analyzed with an irregular wavelet
than with a smooth sinusoid, just as some foods are better handled with a fork
than a spoon. It also makes sense that local features can be described better
with wavelets that have local extent.


2.8 THE CONTINUOUS WAVELET TRANSFORM:


     Mathematically, the process of Fourier analysis is represented by the
Fourier transform:




which is the sum over all time of the signal f(t) multiplied by a complex
exponential. (Recall that a complex exponential can be broken down into real
and imaginary sinusoidal components.) The results of the transform are the
Fourier coefficients F(w), which when multiplied by a sinusoid of frequency w
yields the constituent sinusoidal components of the original signal. Graphically,
the process looks like:




                                  Figure 9
      Similarly, the continuous wavelet transform (CWT) is defined as the sum
over all time of the signal multiplied by scaled, shifted versions of the wavelet




     The result of the CWT is a series many wavelet coefficients C, which are
a function of scale and position.


Multiplying each coefficient by the appropriately scaled and shifted wavelet
yields the constituent wavelets of the original signal:




                                    Figure 10




2.8.1 Scaling:
            We‘ve already alluded to the fact that wavelet analysis produces a
time-scale view of a signal and now we‘re talking about scaling and shifting
wavelets.
What exactly do we mean by scale in this context?
Scaling a wavelet simply means stretching (or compressing) it.


To go beyond colloquial descriptions such as ―stretching,‖ we introduce the
scale factor, often denoted by the letter a.
If we‘re talking about sinusoids, for example the effect of the scale factor is
very easy to see:




                                   Figure 11
The scale factor works exactly the same with wavelets. The smaller the scale
factor, the more ―compressed‖ the wavelet.




                                   Figure 12
It is clear from the diagrams that for a sinusoid sin (wt) the scale factor ‗a‘ is
related (inversely) to the radian frequency ‗w‘. Similarly, with wavelet analysis
the scale is related to the frequency of the signal.
2.8.2 Shifting:


      Shifting a wavelet simply means delaying (or hastening) its onset.
Mathematically, de                              k                            -k)




                                  Figure 13


Five Easy Steps to a Continuous Wavelet Transform:


                     The continuous wavelet transform is the sum over all
time of the signal multiplied by scaled, shifted versions of the wavelet. This
process produces wavelet coefficients that are a function of scale and position.


It‘s really a very simple process. In fact, here are the five steps of an easy
recipe for creating a CWT:


1. Take a wavelet and compare it to a section at the start of the original signal.


2. Calculate a number C that represents how closely correlated the wavelet is
with this                  section of the signal. The higher C is, the more the
similarity. More precisely, if the signal energy and the wavelet energy are
equal to one, C may be interpreted as a      correlation coefficient.
Note that the results will depend on the shape of the wavelet you choose.
                                  Figure 14


3. Shift the wavelet to the right and repeat steps 1 and 2 until you‘ve covered
the whole signal.




                                  Figure 15


4. Scale (stretch) the wavelet and repeat steps 1 through 3.
                                  Figure 16
5. Repeat steps 1 through 4 for all scales.


     When you‘re done, you‘ll have the coefficients produced at different
scales by different sections of the signal. The coefficients constitute the
results of a regression of the original signal performed on the wavelets.


     How to make sense of all these coefficients? You could make a plot on
which the x-axis represents position along the signal (time), the y-axis
represents scale, and the color at each x-y point represents the magnitude of
the wavelet coefficient C. These are the coefficient plots generated by the
graphical tools.




                                  Figure 17


These coefficient plots resemble a bumpy surface viewed from above.
If you could look at the same surface from the side, you might see something
like this:




                                  Figure 18
The continuous wavelet transform coefficient plots are precisely the time-
scale view of the signal we referred to earlier. It is a different view of signal
data than the time- frequency Fourier view, but it is not unrelated.


2.8.3 Scale and Frequency:


      Notice that the scales in the coefficients plot (shown as y-axis labels)
run from 1 to 31. Recall that the higher scales correspond to the most
―stretched‖ wavelets. The more stretched the wavelet, the longer the portion
of the signal with which it is being compared, and thus the coarser the signal
features being measured by the wavelet coefficients.




                                  Figure 19
Thus, there is a correspondence between wavelet scales and frequency as
revealed by wavelet analysis:
• Low scale a=> Compressed wavelet => Rapidly changing details => High
 frequency ‗w‘.
• High scale a=>Stretched wavelet=>Slowly changing, coarse features=>Low
 frequency ‗w‘.
2.8.4 The Scale of Nature:


      It‘s important to understand the fact that wavelet analysis does not
produce a time-frequency view of a signal is not a weakness, but a strength of
the technique.
Not only is time-scale a different way to view data, it is a very natural way to
view data deriving from a great number of natural phenomena.


      Consider a lunar landscape, whose ragged surface (simulated below) is
a result of centuries of bombardment by meteorites whose sizes range from
gigantic boulders to dust specks.
      If we think of this surface in cross-section as a one-dimensional signal,
then it is reasonable to think of the signal as having components of different
scales—large features carved by the impacts of large meteorites, and finer
features abraded by small meteorites.




                                    Figure 20


Here is a case where thinking in terms of scale makes much more sense
than thinking in terms of frequency. Inspection of the CWT coefficients plot for
this signal reveals patterns among scales and shows the signal‘s possibly
fractal nature.
                                  Figure 21
     Even though this signal is artificial, many natural phenomena — from
the intricate branching of blood vessels and trees, to the jagged surfaces of
mountains and fractured metals — lend themselves to an analysis of scale.




2.9 THE DISCRETE WAVELET TRANSFORM:
      Calculating wavelet coefficients at every possible scale is a fair amount
of work, and it generates an awful lot of data. What if we choose only a subset
of scales and positions at which to make our calculations? It turns out rather
remarkably that if we choose scales and positions based on powers of two—
so-called dyadic scales and positions—then our analysis will be much more
efficient and just as accurate. We obtain such an analysis from the discrete
wavelet transform (DWT).


      An efficient way to implement this scheme using filters was developed
in 1988 by Mallat. The Mallat algorithm is in fact a classical scheme known in
the signal processing community as a two-channel sub band coder. This very
practical filtering algorithm yields a fast wavelet transform — a box into which
a signal passes, and out of which wavelet coefficients quickly emerge. Let‘s
examine this in more depth.
2.10 ONE-STAGE FILTERING:
2.10.1 Approximations and Details:


       For many signals, the low-frequency content is the most important part.
It is what gives the signal its identity. The high-frequency content on the other
hand imparts flavor or nuance. Consider the human voice. If you remove the
high-frequency components, the voice sounds different but you can still tell
what‘s being said. However, if you remove enough of the low-frequency
components, you hear gibberish. In wavelet analysis, we often speak of
approximations and details. The approximations are the high-scale, low-
frequency components of the signal. The details are the low-scale, high-
frequency components.
The filtering process at its most basic level looks like this:




                                        Figure 23
       The original signal S passes through two complementary filters and
emerges as two signals.
       Unfortunately, if we actually perform this operation on a real digital
signal, we wind up with twice as much data as we started with. Suppose, for
instance that the original signal S consists of 1000 samples of data. Then the
resulting signals will each have 1000 samples, for a total of 2000.
       These signals A and D are interesting, but we get 2000 values instead
of the 1000 we had. There exists a more subtle way to perform the
decomposition using wavelets. By looking carefully at the computation, we
may keep only one point out of two in each of the two 2000-length samples to
get the complete information. This is the notion of own sampling. We produce
two sequences called cA and cD.




                                      Figure 24


       The process on the right which includes down sampling produces
DWT
Coefficients. To gain a better appreciation of this process let‘s perform a one-
stage discrete wavelet transform of a signal. Our signal will be a pure sinusoid
with
high- frequency noise added to it.


Here is our schematic diagram with real signals inserted into it:




                                  Figure 25
The MATLAB code needed to generate s, cD, and cA is:
s = sin(20*linspace(0,pi,1000)) + 0.5*rand(1,1000);
[cA, cD] = dwt(s,'db2');


       where db2 is the name of the wavelet we want to use for the analysis.
Notice that the detail coefficients cD is small and consist mainly of a high-
frequency noise, while the approximation coefficients cA contains much less
noise than does the original signal.


[length(cA) length(cD)]
ans = 501 501
      You may observe that the actual lengths of the detail and approximation
coefficient vectors are slightly more than half the length of the original signal.
This has to do with the filtering process, which is implemented by convolving
the signal with a filter. The convolution ―smears‖ the signal, introducing
several extra samples into the result.




2.10.2 Multiple-Level Decomposition:


       The decomposition process can be iterated, with successive
approximations being decomposed in turn, so that one signal is broken down
into many lower resolution components. This is called the wavelet
decomposition tree.




                                  Figure 26
Looking at a signal‘s wavelet decomposition tree can yield valuable
information.




                                  Figure 27




2.10.3 Number of Levels:


       Since the analysis process is iterative, in theory it can be continued
indefinitely. In reality, the decomposition can proceed only until the individual
details consist of a single sample or pixel. In practice, you‘ll select a suitable
number of levels based on the nature of the signal, or on a suitable criterion
such as entropy.


2.10.4 Wavelet Reconstruction:


       We‘ve learned how the discrete wavelet transform can be used to
analyze or decompose, signals and images. This process is called
decomposition or analysis. The other half of the story is how those
components can be assembled back into the original signal without loss of
information. This process is called reconstruction, or synthesis. The
mathematical manipulation that effects synthesis is called the inverse discrete
wavelet transforms (IDWT). To synthesize a signal in the Wavelet Toolbox,
we reconstruct it from the wavelet coefficients:




                                   Figure 28


       Where wavelet analysis involves filtering and down sampling, the
wavelet reconstruction process consists of up sampling and filtering. Up
sampling is the process of lengthening a signal component by inserting zeros
between samples:




                                   Figure 29
       The Wavelet Toolbox includes commands like idwt and waverec that
perform   single-level   or   multilevel   reconstruction   respectively on   the
components of one-dimensional signals. These commands have their two-
dimensional analogs, idwt2 and waverec2.
2.10.5 Reconstruction Filters:
      The filtering part of the reconstruction process also bears some
discussion, because it is the choice of filters that is crucial in achieving perfect
reconstruction of the original signal. The down sampling of the signal
components performed during the decomposition phase introduces a
distortion called aliasing. It turns out that by carefully choosing filters for the
decomposition and reconstruction phases that are closely related (but not
identical), we can ―cancel out‖ the effects of aliasing.
      The low- and high pass decomposition filters (L and H), together with
their associated reconstruction filters (L' and H'), form a system of what is
called quadrature mirror filters:




                                    Figure 30


2.10.6 Reconstructing Approximations and Details:


      We have seen that it is possible to reconstruct our original signal from
the coefficients of the approximations and details.




                                       Figure 31
          It is also possible to reconstruct the approximations and details
themselves from their coefficient vectors.
          As an example, let‘s consider how we would reconstruct the first-level
approximation A1 from the coefficient vector cA1. We pass the coefficient
vector cA1 through the same process we used to reconstruct the original
signal. However, instead of combining it with the level-one detail cD1, we feed
in a vector of zeros in place of the detail coefficients
vector:




                                     Figure 32


          The process yields a reconstructed approximation A1, which has the
same length as the original signal S and which is a real approximation of it.
Similarly, we can reconstruct the first-level detail D1, using the analogous
process:




                                     Figure 33
          The reconstructed details and approximations are true constituents of
the original signal. In fact, we find when we combine them that:
                       A1 + D1 = S
       Note that the coefficient vectors cA1 and cD1—because they were
produced by Down sampling and are only half the length of the original signal
— cannot directly be combined to reproduce the signal.
       It is necessary to reconstruct the approximations and details before
combining them. Extending this technique to the components of a multilevel
analysis, we find that similar relationships hold for all the reconstructed signal
constituents. That is, there are several ways to reassemble the original signal:




                                  Figure 34
2.11 RELATIONSHIP OF FILTERS TO WAVELET SHAPES:
       In the section ―Reconstruction Filters‖, we spoke of the importance of
choosing the right filters. In fact, the choice of filters not only determines
whether perfect reconstruction is possible, it also determines the shape of the
wavelet we use to perform the analysis. To construct a wavelet of some
practical utility, you seldom start by drawing a waveform. Instead, it usually
makes more sense to design the appropriate quadrature mirror filters, and
then use them to create the waveform. Let‘s see
how this is done by focusing on an example.
Consider the low pass reconstruction filter (L') for the db2 wavelet.
                          Wavelet function position




                                  Figure 35
The filter coefficients can be obtained from the dbaux command:


Lprime = dbaux(2)


Lprime = 0.3415 0.5915 0.1585 –0.0915


If we reverse the order of this vector (see wrev), and then multiply every even


sample by –1, we obtain the high pass filter H':


Hprime = –0.0915 –0.1585 0.5915 –0.3415


Next, up sample Hprime by two (see dyad up), inserting zeros in alternate


Positions:


HU =–0.0915 0 –0.1585 0 0.5915 0 –0.3415 0


Finally, convolve the up sampled vector with the original low pass filter:


H2 = conv(HU, Lprime);


plot(H2)




                                  Figure 36
       If we iterate this process several more times, repeatedly up sampling
and convolving the resultant vector with the four-element filter vector Lprime,
a pattern begins to emerge:




                                  Figure 37
       The curve begins to look progressively more like the db2 wavelet. This
means that the wavelet‘s shape is determined entirely by the coefficients of
the reconstruction filters. This relationship has profound implications. It means
that you cannot choose just any shape, call it a wavelet, and perform an
analysis. At least, you can‘t choose an arbitrary wavelet waveform if you want
to be able to reconstruct the original signal accurately. You are compelled to
choose a shape determined by quadrature mirror decomposition filters.


2.11.1 The Scaling Function:


     We‘ve seen the interrelation of wavelets and quadrature mirror filters.
The wavelet function       is determined by the high pass filter, which also
produces the details of the wavelet decomposition.
       There is an additional function associated with some, but not all
wavelets. This is the so-called scaling function . The scaling function is very
similar to the wavelet function. It is determined by the low pass quadrature
mirror filters, and thus is associated with the approximations of the wavelet
decomposition. In the same way that iteratively up- sampling and convolving
the high pass filter produces a shape approximating the wavelet function,
iteratively up-sampling and convolving the low pass filter produces a shape
approximating the scaling function.


2.11.2 Multi-step Decomposition and Reconstruction:
A multi step analysis-synthesis process can be represented as:




                                 Figure 38


     This process involves two aspects: breaking up a signal to obtain the
wavelet coefficients, and reassembling the signal from the coefficients. We‘ve
already discussed decomposition and reconstruction at some length. Of
course, there is no point breaking up a signal merely to have the satisfaction
of immediately reconstructing it. We may modify the wavelet coefficients
before performing the reconstruction step. We perform wavelet analysis
because the coefficients thus obtained have many known uses, de-noising
and compression being foremost among them. But wavelet analysis is still a
new and emerging field. No doubt, many uncharted uses of the wavelet
coefficients lie in wait. The Wavelet Toolbox can be a means of exploring
possible uses and hitherto unknown applications of wavelet analysis. Explore
the toolbox functions and see what you discover.
2.11.3 WAVELET DECOMPOSITION:
             Images are treated as two dimensional signals, they change
horizontally and vertically, thus 2D wavelet analysis must be used for images.
2D wavelet analysis uses the same ‘mother wavelets‘ but requires an extra
step at every level of decomposition. The 1D analysis filtered out the high
frequency information from the low frequency information at every level of
decomposition; so only two sub signals were produced at each level.
         In 2D, the images are considered to be matrices with N rows and M
columns. At every level of decomposition the horizontal data is filtered, then
the approximation and details produced from this are filtered on columns.




                  Fig 1: Decomposition of an Image
          At every level, four sub-images are obtained; the approximation, the
vertical detail, the horizontal detail and the diagonal detail. Below the Saturn
image has been decomposed to one level. The wavelet analysis has found
how the image changes vertically, horizontally and diagonally.




       Fig 2:2-D Decomposition of Saturn Image to level 1


To get the next level of decomposition the approximation sub-image is
decomposed, this idea can be seen in figure 3.
           Fig 3: Saturn Image decomposed to Level 3.
Only the 9 detail sub-images and the final sub-image is required to
reconstruct the image perfectly.


When compressing with orthogonal wavelets the energy retained is:
The number of zeros in percentage is defined by:




   Figure: Information flow diagram in image fusion scheme
             employing multi-scale decomposition.

The information flow diagram of wavelet- based image fusion algorithm is

shown in above figure. In wavelet image fusion scheme, the source images I1

(x,y) and I2 (x,y), are decomposed into approximation and detailed

coefficients at required level using DWT. The approximation and detailed

coefficients of both images are combined using fusion rule Φ.




The fused image (If (x, y)) could be obtained by taking the inverse discrete

wavelet transform (IDWT) as:
     If (x, y) = IDWT [Φ{ DWT (I1 (x, y)), DWT (I2 (x, y))}]    ………. (5)



The fusion rule used in this project is simply averages the approximation

coefficients and picks the detailed coefficient in each sub band with the

largest magnitude.



2.11.4 Entropy:

       Entropy of grayscale image



 E = entropy (I)



       Entropy is a statistical measure of randomness that can be used to

characterize the texture of the input image. Entropy is defined as -sum

(p.*log2 (p)) where p contains the histogram counts returned from imhist. By

default, entropy uses two bins for logical arrays and 256 bins for uint8, uint16,

or double arrays. I can be a multidimensional image. If I have more than two

dimensions, the entropy function treats it as a multidimensional grayscale

image and not as an RGB image. Image can be logical, uint8, uint16, or

double and must be real, nonempty, and no sparse. E is double. Entropy

converts any class other than logical to uint8 for the histogram count

calculation so that the pixel values are discrete and directly correspond to a

bin value.
                 3 INTRODUCTION TO MATLAB




3.1 WHAT IS MATLAB?


      MATLAB® is a high-performance language for technical computing. It

integrates computation, visualization, and programming in an easy-to-use

environment where problems and solutions are expressed in familiar

mathematical notation. Typical uses include


   1. Math and computation

   2. Algorithm development

   3. Data acquisition

   4. Modeling, simulation, and prototyping

   5. Data analysis, exploration, and visualization

   6. Scientific and engineering graphics

   7. Application development, including graphical user interface building.


      MATLAB is an interactive system whose basic data element is an array

that does not require dimensioning. This allows you to solve many technical

computing problems, especially those with matrix and vector formulations, in a

fraction of the time it would take to write a program in a scalar non interactive

language such as C or FORTRAN.
The name MATLAB stands for matrix laboratory. MATLAB was originally

written to provide easy access to matrix software developed by the LINPACK

and EISPACK projects. Today, MATLAB engines incorporate the LAPACK

and BLAS libraries, embedding the state of the art in software for matrix

computation.


      MATLAB has evolved over a period of years with input from many users.

In university environments, it is the standard instructional tool for introductory

and advanced courses in mathematics, engineering, and science. In industry,

MATLAB is the tool of choice for high-productivity research, development, and

analysis.


      MATLAB features a family of add-on application-specific solutions

called toolboxes. Very important to most users of MATLAB, toolboxes allow

you to learn and apply specialized technology. Toolboxes are comprehensive

collections of MATLAB functions (M-files) that extend the MATLAB

environment to solve particular classes of problems. Areas in which toolboxes

are available include signal processing, control systems, neural networks,

fuzzy logic, wavelets, simulation, and many others.


3.2 THE MATLAB SYSTEM:


       The MATLAB system consists of five main parts:


3.2.1 Development Environment:


       This is the set of tools and facilities that help you use MATLAB

functions and files. Many of these tools are graphical user interfaces. It
includes the MATLAB desktop and Command Window, a command history,

an editor and debugger, and browsers for viewing help, the workspace, files,

and the search path.


3.2.2 The MATLAB Mathematical Function:


      This is a vast collection of computational algorithms ranging from

elementary functions like sum, sine, cosine, and complex arithmetic, to more

sophisticated functions like matrix inverse, matrix eigen values, Bessel

functions, and fast Fourier transforms.


3.2.3 The MATLAB Language:


      This is a high-level matrix/array language with control flow statements,

functions, data structures, input/output, and object-oriented programming

features. It allows both "programming in the small" to rapidly create quick and

dirty throw-away programs, and "programming in the large" to create

complete large and complex application programs.


3.2.4 Graphics:


       MATLAB has extensive facilities for displaying vectors and matrices as

graphs, as well as annotating and printing these graphs. It includes high-level

functions for two-dimensional and three-dimensional data visualization, image

processing, animation, and presentation graphics. It also includes low-level

functions that allow you to fully customize the appearance of graphics as well

as to build complete graphical user interfaces on your MATLAB applications.
3.2.5 The MATLAB Application Program Interface (API):


       This is a library that allows you to write C and Fortran programs that

interact with MATLAB. It includes facilities for calling routines from MATLAB

(dynamic linking), calling MATLAB as a computational engine, and for reading

and writing MAT-files.


3.3 MATLAB WORKING ENVIRONMENT:


3.3.1 MATLAB Desktop:


      Matlab Desktop is the main Matlab application window. The desktop

contains five sub windows, the command window, the workspace browser, the

current directory window, the command history window, and one or more

figure windows, which are shown only when the user displays a graphic.


      The command window is where the user types MATLAB commands

and expressions at the prompt (>>) and where the output of those commands

is displayed. MATLAB defines the workspace as the set of variables that the

user creates in a work session. The workspace browser shows these

variables and some information about them. Double clicking on a variable in

the workspace browser launches the Array Editor, which can be used to

obtain information and income instances edit certain properties of the variable.


     The current Directory tab above the workspace tab shows the contents

of the current directory, whose path is shown in the current directory window.

For example, in the windows operating system the path might be as follows:

C:\MATLAB\Work, indicating that directory ―work‖ is a subdirectory of the main
directory ―MATLAB‖; WHICH IS INSTALLED IN DRIVE C. clicking on the

arrow in the current directory window shows a list of recently used paths.

Clicking on the button to the right of the window allows the user to change the

current directory.


       MATLAB uses a search path to find M-files and other MATLAB related

files, which are organize in directories in the computer file system. Any file

run in MATLAB must reside in the current directory or in a directory that is on

search path. By default, the files supplied with MATLAB and math works

toolboxes are included in the search path. The easiest way to see which

directories are on the search path. The easiest way to see which directories

are soon the search path, or to add or modify a search path, is to select set

path from the File menu the desktop, and then use the set path dialog box. It

is good practice to add any commonly used directories to the search path to

avoid repeatedly having the change the current directory.


      The Command History Window contains a record of the commands a

user has entered in the command window, including both current and

previous MATLAB sessions. Previously entered MATLAB commands can be

selected and re-executed from the command history window by right clicking

on a command or sequence of commands. This action launches a menu from

which to select various options in addition to executing the commands. This is

useful to select various options in addition to executing the commands. This is

a useful feature when experimenting with various commands in a work

session.
3.3.2 Using the MATLAB Editor to create M-Files:


      The MATLAB editor is both a text editor specialized for creating M-files

and a graphical MATLAB debugger. The editor can appear in a window by

itself, or it can be a sub window in the desktop. M-files are denoted by the

extension .m, as in pixelup.m. The MATLAB editor window has numerous

pull-down menus for tasks such as saving, viewing, and debugging files.

Because it performs some simple checks and also uses color to differentiate

between various elements of code, this text editor is recommended as the tool

of choice for writing and editing M-functions. To open the editor , type edit at

the prompt opens the M-file filename‘s in an editor window, ready for editing.

As noted earlier, the file must be in the current directory, or in a directory in

the search path.


3.3.3 Getting Help:


     The principal way to get help online is to use the MATLAB help browser,

opened as a separate window either by clicking on the question mark symbol

(?) on the desktop toolbar, or by typing help browser at the prompt in the

command window. The help Browser is a web browser integrated into the

MATLAB desktop that displays a Hypertext Markup Language(HTML)

documents. The Help Browser consists of two panes, the help navigator pane,

used to find information, and the display pane, used to view the information.

Self-explanatory tabs other than navigator pane are used to perform a search.
                       4 SOURCE CODE
function varargout = aaaa(varargin)
% AAAA M-file for aaaa.fig
%      AAAA, by itself, creates a new AAAA or raises the existing
%      singleton*.
%
%      H = AAAA returns the handle to a new AAAA or the handle to
%      the existing singleton*.
%
%      AAAA('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in AAAA.M with the given input
arguments.
%
%      AAAA('Property','Value',...) creates a new AAAA or raises the
%      existing singleton*. Starting from the left, property value
pairs are
%      applied to the GUI before aaaa_OpeningFunction gets called.
An
%      unrecognized property name or invalid value makes property
application
%      stop. All inputs are passed to aaaa_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu. Choose "GUI allows
only one
%      instance to run (singleton)".
%
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton', gui_Singleton, ...
                   'gui_OpeningFcn', @aaaa_OpeningFcn, ...
                   'gui_OutputFcn', @aaaa_OutputFcn, ...
                   'gui_LayoutFcn', [] , ...
                   'gui_Callback',   []);
if nargin & isstr(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
     [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
     gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT


% --- Executes just before aaaa is made visible.
function aaaa_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to aaaa (see VARARGIN)

% Choose default command line output for aaaa
handles.output = hObject;
A=ones(256,256);
axes(handles.axes1);
imshow(A);
axes(handles.axes2);

imshow(A);
axes(handles.axes3);
imshow(A);
% Update handles structure
guidata(hObject, handles);

% UIWAIT makes aaaa wait for user response (see UIRESUME)
% uiwait(handles.figure1);


% --- Outputs from this function are returned to the command line.
function varargout = aaaa_OutputFcn(hObject, eventdata, handles)
% varargout cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
varargout{1} = handles.output;


% --- Executes on button press in image1.
function image1_Callback(hObject, eventdata, handles)
% hObject    handle to image1 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

[filename, pathname] = uigetfile('*.jpg', 'Pick an Image');

       if isequal(filename,0) | isequal(pathname,0)

          warndlg('User pressed cancel')
       else

M1 = imread( filename);
       handles.M1 = M1;
       axes(handles.axes1);
       imshow(M1);

 end

[r c p]=size(M1);
if p==3
    M1=rgb2gray(M1);
end
M1=imresize(M1,[256,256]);
%
        %Update handles structure
            handles.M1=M1;
        guidata(hObject, handles);
% --- Executes on button press in image2.
function image2_Callback(hObject, eventdata, handles)
% hObject    handle to image2 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

[filename, pathname] = uigetfile('*.jpg', 'Pick an Image');

    if isequal(filename,0) | isequal(pathname,0)

       warndlg('User pressed cancel')
    else

       [M2,ma] = imread([pathname, filename]);
       handles.M2 = M2;
       axes(handles.axes2);
       imshow(M2);
     end
     [r c p]=size(M2);
     if p==3
         M2=rgb2gray(M2);
     end
     M2=imresize(M2,[256,256]);

       handles.M2=M2;
       %Update handles structure
       guidata(hObject, handles);


% --- Executes on button press in image_fuse.
function image_fuse_Callback(hObject, eventdata, handles)
% hObject    handle to image_fuse (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

filename2 = handles.M1;
filename1 = handles.M2;
contents = get(handles.popupmenu1,'Value');

switch contents

    case 1
        M1=filename1;
        M2=filename2;
        M1=double(M1);
        M2=double(M2);
        F = selb(M1,M2);
        axes(handles.axes3);
        imshow(F,[]);
        title('Averaging');

        e = entropy(F);
        disp('entrophy');
        disp(e);

    case 2
        M1=filename1;
        M2=filename2;
        M1=double(M1);
          M2=double(M2);
          F = selc(M1,M2);
          axes(handles.axes3);
          imshow(F,[]);

          e = entropy(F);
          disp('entrophy');
          disp(e);

      case 3
          M1=filename1;
          M2=filename2;
          M1=double(M1);
          M2=double(M2);
          F = -selc(-M1,-M2);
          axes(handles.axes3);
          imshow(F,[]);

          e = entropy(F);
          disp('entrophy');
          disp(e);

      case 4
          M1=filename1;
          M2=filename2;
          M1=double(M1);
          M2=double(M2);
          Y = fuse_pca(M1,M2)
          axes(handles.axes3);
          imshow(Y,[]);

          e = entropy(Y);
          disp('entrophy');
          disp(e);
end

% --- Executes during object creation, after setting all properties.
function popupmenu1_CreateFcn(hObject, eventdata, handles)
% hObject    handle to popupmenu1 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns
called

% Hint: popupmenu controls usually have a white background on
Windows.
%        See ISPC and COMPUTER.
if ispc
     set(hObject,'BackgroundColor','white');
else

set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'
));
end


% --- Executes on selection change in popupmenu1.
function popupmenu1_Callback(hObject, eventdata, handles)
% hObject    handle to popupmenu1 (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles     structure with handles and user data (see GUIDATA)

% Hints: contents = get(hObject,'String') returns popupmenu1 contents
as cell array
%        contents{get(hObject,'Value')} returns selected item from
popupmenu1


% --- Executes on button press in BACK.
function BACK_Callback(hObject, eventdata, handles)
% hObject    handle to BACK (see GCBO)
% eventdata reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

close aaaa;
% gui_main;




                                   \
              5 OUTPUTS




             5.1 Fig: Image 1




             5.2 Fig: Image 2




5.3 Figure: Fused image by simple average.
  5.4 Figure: Fused image by PCA




5.5 Figure: Fused image by Wavelets
                          6 CONCLUSIONS




      Pixel-level image fusion using wavelet transform and principal

component analysis are implemented in PC MATLAB. Different image fusion

performance metrics with and without reference image have been evaluated.

The simple averaging fusion algorithm shows degraded performance. Image

fusion using wavelets with higher level of decomposition shows better

performance in some metrics while in other metrics, the PCA shows better

performance. Some further investigation is needed to resolve this issue.
                               7 REFERENCES



1. Gonzalo, Pajares & Jesus Manuel, de la Cruz. A wavelet-based image

fusion tutorial. Pattern Recognition, 2004, 37, 1855-872.

2. Varsheny, P.K. Multi-sensor data fusion. Elec. Comm. Engg., 1997, 9(12),

245-53.

3. Mallet, S.G. A theory for multiresolution signal decomposition: The wavelet

representation. IEEE Trans. Pattern Anal. Mach. Intel., 1989, 11(7), 674-93.

4. Wang, H.; Peng, J. & Wu, W. Fusion algorithm for multisensor image based

on discrete multi wavelet transform. IEE Proc. Visual Image Signal Process.,

2002, 149(5).

5. Mitra Jalili-Moghaddam. Real-time multi-focus image fusion using discrete

wavelet transform and Laplasican pyramid transform. Chalmess University of

Technology, Goteborg, Sweden, 2005. Masters thesis.

6. Daubechies, I. Ten lectures on wavelets. In Regular Conference Series in

Applied Math‘s, Vol. 91, 1992, SIAM, Philadelphia.

7. h t t p : / / e n . w i k i p e d i a . o r g / w i k i / Principal components analysis.

8. Naidu, V.P.S.; Girija, G. & Raol, J.R. Evaluation of data association and

fusion algorithms for tracking in the presence of measurement loss. In AIAA

Conference on Navigation, Guidance and Control, Austin, USA, August 2003,

pp. 11-14.

9. Arce, Gonzalo R. Nonlinear signal processing A statistical approach. Wiley-

Interscience Inc. Publication, USA, 2005.

9. Arce, Gonzalo R. Nonlinear signal processing A statistical approach. Wiley-

Interscience Inc. Publication, USA, 2005.

				
DOCUMENT INFO