\usepackage{.\msuthesis,psfig}

\documentclass[11pt]{report}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\begin{document}

\title{INTEGRATING MODELING AND MEASUREMENT IN THE MAPPING OF SIGNAL PROCESSING
ALGORITHMS TO ADVANCED ARCHITECTURES}
\author{Nayda G. Santiago Santiago}
\date{June 1999}
\maketitle

\begin{abstract}
Signal processing applications are frequently computationally intensive

tasks. When developing algorithms for these applications, implementation

issues should be considered in order to obtain efficient coding. However,

this imposes an additional burden on the signal processing specialist who,

in most cases, is more concerned about correctness and effectiveness issues.

This proposal discusses studying the mapping process of algorithms to high

performance platforms, to correlate the different options that might be

available to the signal processing specialist. As an outcome, it is expected

to produce high-level representations of the available mapping space,

leading to the selection of the best alternatives to obtain an efficient

implementation. To effectively carry out this study, it is proposed to

incorporate both system modeling and measurement in the process of

implementing algorithms on different advanced architectures. A plan for the

review of related work, problem statement, plan for the research, and

timeline with the expected duration of the different activities are

presented.
\end{abstract}

\tableofcontents

\article{A DISSERTATION PROPOSAL}

\school{Michigan State University}

\degree{DOCTOR OF PHILOSOPHY}

\department{Department of Electrical and Computer Engineering}

\chapter{Introduction}

Signal processing research has led to a number of revolutionary applications

which have changed communications, medicine, defense, marketing, and even

entertainment. As a consequence of this trend, media processing has steadily

become one of the dominant computing workloads nowadays \cite{RaAd99}. The

development of algorithms for signal processing applications frequently

involves a process of mathematical study, proof, simulation, and

implementation. Researchers involved in algorithm development are, in most

cases, signal processing specialists who typically have a greater expertise

in the mathematics of the problem than in computing. When a signal

processing specialist makes implementational decisions to select a specific

factor such as a platform, or a computer language, or diverse compiler

options, or a programming paradigm, among others, the performance of the

implementation is affected.

A common process followed in the implementation of an algorithm is going

from the highest level of abstraction to mapping the procedure into the

selected target architecture. If the algorithm is inherently complex, it

could require implementation on a high performance system, extending the

mapping process. Further, new algorithm implementations should have improved

qualities, such as efficiency and performance. All these requirements, when

combined together, make the mapping process a challenging task, calling for

the need of friendly tools to assist the signal processing specialist.

Problem solving environments (PSEs) are an emerging area in computational

science and engineering aimed at reducing the turnaround time in the

production of efficient code for computationally intensive scientific

applications [2-8].\nocite{GaRi94} \nocite{RiBo96} \nocite{MaWe97} \nocite
{applications [2-8].\nocite{GaRi94} \nocite{RiBo96} \nocite{MaWe97} \nocite}{%
AkKa97} \nocite{EnWo97} \nocite{FoIl87} \nocite{Felm87} According to Rice

and Boisvert, a problem solving environment is an integrated system that

contains all the facilities required to solve a specific class of problems

in an efficient and convenient way \cite{RiBo96}. The availability of such

type of systems for the type of problems found in the area of signal

processing would provide a great aid in the development and implementation

of new algorithms. Several attempts have been reported towards the

development of such type of tools [9-11].\nocite{ShMo94} \nocite{BoFr97}

\nocite{PiHa95} However, none of them have effectively addresses the problem

of producing a friendly, efficient tool integrating performance measures and

supporting the development of new algorithms.

The central objective of this research is to study the mapping process of a

class of signal processing applications into high-performance architectures,

and to develop a set of high-level abstractions such that, given a specific

problem and a target architecture, a set of guidelines are produced which

would lead to an efficient implementation. The abstractions resulting from

this research will combine information from performance modeling,

performance monitoring, and mapping techniques in a coherent way such that

they can be either incorporated with an existing PSE or be used as the core

for the development of a new PSE in the selected signal processing area.

The tools to be developed will have as target application the areas of image

processing and pattern recognition. A case study has been selected for the

development of new algorithms for wavelet-based feature extraction and

texture analysis for database image retrieval. These extracted features are

intended to be mapped to the compressed domain and produce faster database

retrieval mechanisms.

\emph{In summary, this research will explore the solution space for the
process of designing and mapping texture analysis algorithms to advanced
architectures. The process is supported by a problem-solving environment,
which increases productivity, and the objective is to improve the efficiency
of the implementation. The solution space consists of platform configuration
options at different levels in the mapping process, performance monitoring
results, and algorithm alternatives. The goal is to determine correlation
and provide appropriate feedback to the signal processing specialist to
guide him/her in algorithmic decisions.}

The rest of this document is organized as follows. The next chapter presents

previous research reported in the areas of this proposal. Chapter 3 contains

a detailed statement of the problem, followed by a research plan and a

timeline in chapter 4. The last chapter contains a discussion of the

expected contributions of this research.

\chapter{Problem Background}

Diverse solutions to the problem of implementing algorithms to different

architectures have been reported. Different environments and tools have been

developed to assist computer programmers in the mapping process. This

chapter presents background information on the previous work reported in

these areas and the signal processing application area addressed by this

research.

To properly understand the mapping process, it should be tied to a specific

application area. The first section of this chapter provides a discussion on

the application addressed in this research: compressed domain techniques for

image processing. It first discusses the current standard JPEG compression

and some compressed domain techniques developed for JPEG. Wavelet-based

compression and texture analysis are discussed next.

A discussion on implementation of algorithms and current advanced

architectures follows. This information provides the background for the

architectures to be used in this research. Finally, a section on

environments and tools related to the proposed work is presented. A

description of the current state of problem solving environments for signal

processing applications, performance models and tools, and the message

passing interface for parallel programming is given.

\section{Compressed Domain Techniques for Image Processing}

Multimedia applications include images and video as standard components.

Most visual data is stored in compressed format to improve storage capacity

and transmission \cite{ShSe95}. Image processing applications requiring

real-time response have two main problems: the volume of data to manipulate

is large and the compression and decompression algorithms are complex. This

is computationally demanding and therefore media processing is becoming one

of the dominant computing workloads of the future \cite{RaAd99}.

Research is currently being developed in the area of image manipulation in

the compressed domain [12-16]. \nocite{ShSe95} \nocite{Ch95a} \nocite{Ch95b}

\nocite{SmRo93} \nocite{Ye96} Compressed domain methods imply that the image

is compressed or only partially decompressed when it is manipulated.

The advantages of compressed domain techniques are mainly three:

\begin{itemize}
\item  Reduced computational cost

\item  Elimination of the overhead of coding/decoding

\item  Some compression algorithms aid in the image content analysis
\end{itemize}

Compressed domain techniques presented in literature are applied to images

compressed in the JPEG (Joint Photographic Experts Group) standard format.

JPEG is based on the computation of the discrete cosine transform. Different

compressed domain techniques can be applied to other compression methods

such as wavelet image compression.

One of the problems faced when retrieving images in databases is the

consistency in the keywords used to describe the image \cite{ChSm95}. These

will depend on the person describing it. Image features such as shape,

color, or texture can be used for retrieving images from databases to solve

this problem.

Image compression using the discrete wavelet transform has been very

effective, especially at low bit rates [18-22]. \nocite{Shap93} \nocite
{effective, especially at low bit rates [18-22]. \nocite{Shap93} \nocite}{%
BuSi97} \nocite{NgWe97} \nocite{RiVe91} \nocite{BuGo98} The use of the

discrete wavelet transform in the analysis of texture has proven successful

in current research work [23-28]. \nocite{ScLi97} \nocite{ChBa98} \nocite{in
current research work [23-28]. \nocite{ScLi97} \nocite{ChBa98} \nocite}{%
PoCa97} \nocite{SaLi95} \nocite{PiTe96} \nocite{VaSc99} The synergy between

image compression and feature extraction using the discrete wavelet

transform can be exploited to improve performance of algorithms \cite{Ch95b}.

The following sections will provide background on current compressed domain

techniques for JPEG and wavelet compressed images.

\subsection{JPEG and compressed domain techniques}

The algorithm for JPEG\ standard image compression starts by partitioning

the image in blocks \cite{SmRo93,Wa91}. Each block is an $8\times 8$ matrix

of values in the range $[0,$ $255]$. The block is normalized to bring the

values to the range $[-128$, $127]$. The Discrete Cosine Transform (DCT) is

computed on each $8\times 8$ matrix using equation (\ref{e:cos1}). Let $

y[i,j]$ denote the input signal. The DCT is defined as

\begin{equation}
Y\left[ u,v\right] =\frac{1}{4}\sum\limits_{i}\sum\limits_{j}C\left(
i,u\right) C\left( j,v\right) y[i,j]  \label{e:cos1}
\end{equation}

where

\begin{equation}
C\left( i,u\right) =A(u)\cos \frac{\left( 2i+1\right) u\pi }{16}
\label{e:cos2}
\end{equation}

\begin{equation}
A(u)=\left\{ 
\begin{array}{cc}
\frac{1}{\sqrt{2}} & u=0 \\ 
1 & \text{otherwise}
\end{array}
\right.
\end{equation}

The lowest frequency components of the DCT are in the upper left corner of

the matrix and the highest frequency components are in the lower right

corner. The next step in JPEG is the quantization of the DCT elements based

on a quantization table. The quality of the image depends on the values in

the quantization table. After the quantization step, the $8\times 8$ matrix

is converted into a $64$-element vector by zig-zag ordering. This ordering

places the lowest frequency component at the beginning of the vector, and

sorts the components in increasing frequency content.

Once the zig-zag ordering is done, the resulting vectors contain a large

number of zeros. Run length encoding (RLE) is applied to encode the vector

into $(skip,value)$ pairs. The final step is to use entropy coding (for

example, Huffman coding) to compress the RLE blocks.

Compressed domain techniques are based on the properties of the JPEG and the

DCT. The properties are studied to map image manipulation from the

space-time domain to the compressed domain. Most of these methods operate in

the transform domain while others operate on the run-length encoded data

(see figure \ref{f:jpegcompr}).

\FRAME{ftbpFU}{423.9375pt}{187.375pt}{0pt}{\Qcb{Relationship between the
transform domain, the RLE\ domain, and the compression algorithm.}}{\Qlb{
f:jpegcompr}}{Figure }{\special{language "Scientific Word";type
"GRAPHIC";maintain-aspect-ratio TRUE;display "PICT";valid_file "T";width
423.9375pt;height 187.375pt;depth 0pt;original-width
563.25pt;original-height 247.75pt;cropleft "0";croptop "1";cropright
"0.9996";cropbottom "0";tempfilename
'z:/propOSAL/Compress.wmf';tempfile-properties "NP";}}

Two examples of compressed domain techniques are the work of Smith and Rowe

\cite{SmRo93} and the work of Shen and Sethi \cite{ShSe96}. Smith and Rowe

present the implementation of four basic algebraic operations in the

compressed domain in \cite{SmRo93}: scalar multiplication, scalar addition,

pixel addition, and pixel multiplication. These operations are applied on

the RLE domain. These basic operations will build the foundation for other

image manipulation methods in the compressed domain.

Shen and Sethi study direct feature extraction from compressed image for

edge detection and coarse classification in \cite{ShSe96}. These are

obtained in the DCT transform domain by developing a model of an edge. The

first feature to compute in order to do classification is the location of

the area of interest. This will provide the blocks containing edges in the

image. Once the blocks are extracted, three parameters are of interest to

detect edges: edge orientation $(\theta )$, edge strength $(h)$, and edge

offset from the center $(d)$. These will determine the edge strength. The

stronger areas are classified as edges by thresholding. Experiments show

that this approach is much faster but not as precise as the spatial domain

edge detection algorithms. This method is good for coarse image

classification.

\subsection{Wavelet-based Compression and Texture Analysis}

Wavelets are functions that are used to hierarchically decompose data or

other functions. These functions satisfy certain mathematical requirements.

This is similar to Fourier decomposition of functions in terms of sine and

cosine functions. However, unlike Fourier analysis, wavelets analyze data

according to scale or resolution. If we look at a signal at a coarse scale,

we can see gross features, while if we look at it at a fine scale, we will

see its details \cite{Gra95}. This characteristic makes wavelets useful for

approximating data with discontinuities like images. The more adapted the

wavelet is to the data, the more sparse the representation in the wavelet

domain is. This makes wavelets well-suited for compression and feature

extraction in images \cite{Gra95}.

In order to explain the wavelet transform of an image, we will start

explaining the wavelet transform of a one-dimensional signal and then expand

to a two-dimensional signal. It should be noted that we are limiting this

work to wavelets in $\mathcal{Z}_{N}$ (space of cyclic signals).

\subsubsection{Wavelet transform of a 1-D signal}

Wavelets are basis functions in continuous time that have the special

feature that all functions are constructed from a single mother wavelet.

This mother wavelet is a small wave (short signal or localized in time). In

the wavelet transform, a function is represented as a superposition of

dilated and translated wavelets. The discrete wavelet transform on $\mathcal{%
\ Z}_{N}$ is a linear operator that transforms a vector into another of the

same length.

There are infinite number of functions that satisfy the properties of

wavelets. They differ in how localized they are in time or frequency. This

means that there are wavelets that emphasize sharp changes of the signal,

other emphasize smooth changes, ones that give more information about

changes in time, and others that give more information about frequency

content of the signal. Daubechies compact support wavelet is very well

known. Compact support implies that it does not extend to $\pm \infty ,$ and

the filter to compute the wavelet is a finite impulse response filter. A

finite impulse response filter (FIR) is a filter whose output is the

convolution of a \textbf{finite} number of filter coefficients with the

input signal.

The computation of the wavelet coefficients consists of two steps: filtering

the signal and decimating the output of the filter. In \emph{decimation by
two}, every other element of the signal is obtained at the output. Filtering

and decimation are repeated recursively to obtain different resolution

levels \cite{VeKo95}. The scheme to compute the wavelet transform is shown

in figure $\ref{f:1d1lbl}$. First, the signal is filtered with $h_{l}$ (low

pass filter) and $h_{h}$ (high pass filter) and then decimated. The wavelet

transform is $w=[w_{l}~w_{h}]$.

\FRAME{ftbpFU}{270.5625pt}{97.8125pt}{0pt}{\Qcb{Block diagram of 1-D wavelet
transform. One level of transform.}}{\Qlb{f:1d1lbl}}{Figure }{\special%
{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio
TRUE;display "USEDEF";valid_file "T";width 270.5625pt;height 97.8125pt;depth
0pt;original-width 268.875pt;original-height 96.1875pt;cropleft "0";croptop
"1";cropright "1";cropbottom "0";tempfilename
'z:/propOSAL/dwt1d1l.wmf';tempfile-properties "NP";}}

So far we have only computed the finest scale of wavelet coefficient. To

compute a coarser scale, we have to repeat the same process for the low

frequency component part of the signal.

\FRAME{ftbpFU}{431.6875pt}{129.5pt}{0pt}{\Qcb{Two level decomposition
wavelet transform.}}{\Qlb{f:lev2}}{Figure }{\special{language "Scientific
Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "PICT";valid_file
"T";width 431.6875pt;height 129.5pt;depth 0pt;original-width
307.0625pt;original-height 91.25pt;cropleft "0";croptop "1";cropright
"1";cropbottom "0";tempfilename
'z:/propOSAL/dwt1d2l.wmf';tempfile-properties "NP";}}

Figure \ref{f:lev2} shows the recursive computation of the wavelet

transform. This process is repeated on the low frequency component of the

signal until a final scale is reached. This is called a \emph{pyramid
decomposition}.

\subsubsection{Wavelet transform of a 2-D signal}

The standard wavelet decomposition of a 2-D signal or image first computes

the one dimensional wavelet transform in each row of the image obtaining two

images. Then the one dimensional wavelet transform of the columns of these

two images is obtained \cite{StDe95a}. This procedure is repeated on the low

frequency image to generate the next coarser scale of wavelet coefficients.

Figure \ref{f:tdwt} illustrates a block diagram showing the computation of

the two-dimensional discrete wavelet transform.

\begin{figure}[tbph]
\newcounter{cms} \setlength{\unitlength}{1mm}
\par
\begin{picture}(120,75)

\put(20,37){\line(1,0){10}}
\put(5,37){\makebox(0,0)[l]{$image$}}

\put(30,35){\line(0,1){20}}
\put(30,55){\vector(1,0){5}}
\put(30,20){\line(0,1){20}}
\put(30,20){\vector(1,0){5}}

\put(35,37){\makebox(0,0)[l]{\emph{filtering rows}}}

\put(35,50){\framebox(10,10){$h_{l}$}}
\put(45,55){\vector(1,0){5}}
\put(50,50){\framebox(10,10){$\downarrow$ 2}}
\put(60,55){\line(1,0){10}}
\put(35,15){\framebox(10,10){$h_{h}$}}
\put(45,20){\vector(1,0){5}}
\put(50,15){\framebox(10,10){$\downarrow$ 2}}
\put(60,20){\line(1,0){10}}

\put(70,55){\line(0,1){10}}
\put(70,45){\line(0,1){10}}
\put(70,65){\vector(1,0){5}}
\put(70,45){\vector(1,0){5}}

\put(75,55){\makebox(0,0)[l]{\emph{filtering columns}}}
\put(75,20){\makebox(0,0)[l]{\emph{filtering columns}}}

\put(75,60){\framebox(10,10){$h_{l}$}}
\put(85,65){\vector(1,0){5}}
\put(75,40){\framebox(10,10){$h_{h}$}}
\put(85,45){\vector(1,0){5}}

\put(90,60){\framebox(10,10){$\downarrow$ 2}}
\put(100,65){\vector(1,0){5}}
\put(90,40){\framebox(10,10){$\downarrow$ 2}}
\put(100,45){\vector(1,0){5}}

\put(108,65){\makebox(0,0)[l]{$I_{ll}$}}
\put(108,45){\makebox(0,0)[l]{$I_{lh}$}}

\put(70,20){\line(0,1){10}}
\put(70,10){\line(0,1){10}}
\put(70,30){\vector(1,0){5}}
\put(70,10){\vector(1,0){5}}

\put(75,25){\framebox(10,10){$h_{l}$}}
\put(85,30){\vector(1,0){5}}
\put(75,5){\framebox(10,10){$h_{h}$}}
\put(85,10){\vector(1,0){5}}

\put(90,25){\framebox(10,10){$\downarrow$ 2}}
\put(100,30){\vector(1,0){5}}
\put(90,5){\framebox(10,10){$\downarrow$ 2}}
\put(100,10){\vector(1,0){5}}

\put(108,30){\makebox(0,0)[l]{$I_{hl}$}}
\put(108,10){\makebox(0,0)[l]{$I_{hh}$}}

\end{picture}
\caption{Block diagram of one level of decomposition of a two-dimensional
wavelet transform.}
\label{f:tdwt}
\end{figure}

The filtering-downsampling steps are repeated on the \emph{ll} coefficients.

This can be repeated $k=\log _2M$ times where $M$ is $\min

\{Num\_rows,Num\_cols\}$ to give a level-$k$ wavelet decomposition.

\subsubsection{Inverse Wavelet Transform}

To obtain the inverse wavelet transform, the columns of the image are

interpolated or upsampled and then filtered with a synthesis filter. Then

the rows are interpolated and filtered. In interpolation by two, a zero is

added after every element of the vector. Figure \ref{f:itdwt} shows the

block diagram of the synthesis stage of the wavelet transform.

\begin{figure}[tbph]
\setlength{\unitlength}{1mm}
\par
\begin{picture}(120,75)

\put(15,65){\makebox(0,0)[l]{$I_{ll}$}}
\put(15,45){\makebox(0,0)[l]{$I_{lh}$}}
\put(15,30){\makebox(0,0)[l]{$I_{lh}$}}
\put(15,10){\makebox(0,0)[l]{$I_{lh}$}}

\put(20,65){\vector(1,0){5}}
\put(20,45){\vector(1,0){5}}

\put(25,55){\makebox(0,0)[l]{\emph{filtering columns}}}
\put(25,20){\makebox(0,0)[l]{\emph{filtering columns}}}

\put(25,60){\framebox(10,10){$\uparrow$ 2}}
\put(35,65){\vector(1,0){5}}
\put(25,40){\framebox(10,10){$\uparrow$ 2}}
\put(35,45){\vector(1,0){5}}

\put(40,60){\framebox(10,10){$g_{l}$}}
\put(50,65){\line(1,0){10}}
\put(40,40){\framebox(10,10){$g_{h}$}}
\put(50,45){\line(1,0){10}}

\put(20,30){\vector(1,0){5}}
\put(20,10){\vector(1,0){5}}

\put(25,25){\framebox(10,10){$\uparrow$ 2}}
\put(35,30){\vector(1,0){5}}
\put(25,5){\framebox(10,10){$\uparrow$ 2}}
\put(35,10){\vector(1,0){5}}

\put(40,25){\framebox(10,10){$g_{l}$}}
\put(50,30){\line(1,0){10}}
\put(40,5){\framebox(10,10){$g_{h}$}}
\put(50,10){\line(1,0){10}}

\put(60,30){\vector(0,-1){7}}
\put(60,10){\vector(0,1){7}}

\put(60,65){\vector(0,-1){6}}
\put(60,45){\vector(0,1){6}}

\put(60,20){\circle{7}}
\put(60,55){\circle{7}}
\put(60,20){\makebox(0,0)[c]{+}}
\put(60,55){\makebox(0,0)[c]{+}}

\put(64,55){\vector(1,0){6}}
\put(64,20){\vector(1,0){6}}

\put(70,50){\framebox(10,10){$\uparrow$ 2}}
\put(80,55){\vector(1,0){5}}
\put(85,50){\framebox(10,10){$g_{l}$}}
\put(95,55){\line(1,0){10}}
\put(70,15){\framebox(10,10){$\uparrow$ 2}}
\put(80,20){\vector(1,0){5}}
\put(85,15){\framebox(10,10){$g_{h}$}}
\put(95,20){\line(1,0){10}}

\put(70,37){\makebox(0,0)[l]{\emph{filtering rows}}}

\put(105,55){\vector(0,-1){14}}
\put(105,20){\vector(0,1){13}}

\put(105,37){\circle{7}}
\put(105,37){\makebox(0,0)[c]{+}}

\put(109,37){\vector(1,0){9}}
\put(119,37){\makebox(0,0)[l]{\emph{image}}}

\end{picture}
\caption{Block diagram of one-level synthesis of a two-dimensional wavelet
transform.}
\label{f:itdwt}
\end{figure}

\subsubsection{Compression}

Some properties of the discrete wavelet transform make it effective for

image compression. First, since the wavelet representation is localized in

both time and frequency, it provides a compact representation of the image.

Wavelet bases can represent singularities so edges are well represented in

the wavelet transform. Second, it clusters insignificant coefficients \cite
{the wavelet transform. Second, it clusters insignificant coefficients \cite}%
{BuGo98}. This can be exploited in the compression scheme like, for example,

the scheme used in embedded zerotree wavelet (EZW) compression \cite{Shap93}

. Third, multiresolution analysis allows fast coarse image browsing.

Compression methods based on wavelets usually outperform DCT-based JPEG,

specially at low bit rates \cite{BuGo98}.

\subsubsection{Texture analysis}

As the size of image databases increases, the need for faster algorithms for

image retrieval become important. The use of keyword description of images

is discouraged by the lack of consistency in these descriptions [34-35].

\nocite{MeKi97} \nocite{SmCh94} A statistical approach using feature

selection and similarity measures is suggested for this purpose \cite{MeKi97}

. Texture is one of the features that distinguish images.

There are three different statistical methods for texture classification

problems: image-based, model-based, and transform-based statistical

algorithms \cite{NaSm93}. This classification is based on how the features

are computed. The purpose of transform-based methods is to determine the

spatial-frequency component of the signal which indicates texture.

Multiresolution analysis is suitable for this purpose and has been used for

texture analysis [23-28].\nocite{ScLi97} \nocite{ChBa98} \nocite{PoCa97}

\nocite{SaLi95} \nocite{PiTe96} \nocite{VaSc99}

Once the image processing algorithm is developed, the next step is

implementation on a specific architecture.

\section{Implementation of Algorithms}

An algorithm is a set of steps to solve a specific problem. It is an

abstraction of a computer program. This abstraction removes the constraints

of a particular hardware or high-level language from the description on how

to solve a given problem \cite{Naps92}.

The use of an algorithm as a high-level of abstraction allows portability

and faster development time of software. With the rapid changes in

computers, algorithms remain static while architectures change. However,

there is a trade-off between high level of abstraction and execution time.

Tuning a program to a particular computer often requires knowledge of the

computer architecture and instruction set architecture.

An implementation is the task of turning an algorithm into a computer

program \cite{Baas78}. There is no standard procedure on how to do this.

There are a set of decisions to make in the process such as which computer,

language, compiler, etc. It could be easy or difficult based on these

decisions. One abstraction can have different implementations.

Once the implementation is done, the first step is to determine if results

are correct \cite{Wick94}. Results should be what was expected. The second

step is to determine if the implementation is efficient. Efficiency is

typically measured in terms of execution time and memory space \cite{Naps92}

. Execution time is mainly determined by the amount of data and number of

operations to perform. A rough estimate of efficiency can therefore be

obtained prior to implementation by doing a complexity analysis, that is,

doing a mathematical analysis of the estimated number of operations given

the amount of data.

The efficiency of an algorithm can also be determined after implementation

using monitoring tools. These will provide detailed information on the

behavior of the algorithm and quantify successes or failures in the

implementation \cite{Dowd93}. Unfortunately, to understand performance

information provided these tools it is required to have knowledge about the

computer architecture . This is in contrast to the concept of abstraction.

Poor memory access patterns or improper compiler optimizations can prevent a

program from being efficient. Performance tools that can relate high-level

abstractions to efficiency are needed.

A methodology for mapping mathematical algorithms to platforms was

introduced in \cite{RoSa93}. This methodology presents the concept of a

computational mathematical environment which basically consists of an input

data set, an output data set, a set of operators acting on this set, and a

set of rules for these operators. In this environment the parameters of the

computing platform are presented to a computational mathematical framework

to assist in the mapping process. We will extend this work to incorporate

performance metrics to the computational mathematical framework.

\section{Advanced Architectures}

The computing power requirements for scientific applications have led to the

development of different approaches for meeting the demands of processing

speed, memory size and speed, and data input/output rates. Increases in

performance have come from several advances \cite{MeCu98}. Some of these

advances include:

\begin{itemize}
\item  Microprocessors have become faster through the use of

instruction-level parallelism, multilevel caches, and faster clock speeds.

\item  Different schemes have been developed for effective interconnection

between processors and memory.

\item  Users have learned how to use multiple processors and deep memory

hierarchies.

\item  Software tools have been improved.
\end{itemize}

Advances in microprocessor technology include the development of the Intel

Pentium family, the Motorola PowerPC family, and the Digital Equipment

Corporation's Alpha 21164 microprocessor \cite{Rafi99}. These are 32 or

64-bits microprocessors with clock speeds higher than 233 MHz.

Taking as an example the Intel Pentium family, it includes the Intel

Pentium, Pentium Pro, and Pentium II. The Pentium is a superscalar level 2

processor. This is, it can execute two instructions per clock cycle. It uses

two 5-stage pipelines, two 8K internal caches, one for code and one for

data, and two translation lookaside buffers (TLB) to handle the pipelines.

It prefetches up to 32 bytes of instructions for branch-prediction \cite{It
prefetches up to 32 bytes of instructions for branch-prediction \cite}{Rafi99%
}. The Pentium Pro is a superscalar level 3 processor. It has three

12-stage pipelines to execute three instructions per clock cycle. This

processor also combines branch prediction, speculative execution, and

analysis for data dependencies to improve program execution over the Pentium

and to prevent stalling. This last feature is called \emph{dynamic execution}

. Finally, the Pentium II combines the characteristics of the Pentium Pro

and the MMX\ (\textbf{m}atrix \textbf{m}ath e\textbf{x}tensions) technology.

MMX technology is intended for accelerating multimedia and communications

applications \cite{PeWi97}. It contains special instructions aimed at

efficiently manipulating video, audio, and graphics by exploiting the

inherent parallelism present in these applications. These instructions

operate on the data in a single-instruction multiple-data fashion. Like the

Pentium Pro, the Pentium II uses dynamic execution and dual independent bus

architecture (DIB). DIB consists of two independent buses, one for cache and

the other for the system to accelerate the access of information.

Scalable parallel processors (SPP) are also used to meet performance

demands. SPPs are hundreds or thousands of state-of-the-art interconnected

processors. The goal of SPPs is to obtain fast computers through highly

parallel designs and substantial parallelism within the processor. Issues to

consider when designing SPPs are whether to use Single-Instruction

Multiple-Data (SIMD), Multiple-Instruction Multiple-Data (MIMD), and Very

Long Instruction Word (VLIW). In the SIMD architecture, the same instruction

is executed on multiple processors at the same time. In MIMD, different

instructions and data are operated on by each processor. In the VLIW

architecture, multiple functional units are exercised by the instructions on

each cycle. The instruction word is composed of simpler instructions which

are executed simultaneously. The real challenge in SPPs is how to perform at

the potential the system has. Some of the challenges include decisions on

the interconnection hardware, memory system design, compilers, and algorithm

design \cite{MeCu98}.

Interconnection hardware could be a custom design, a commercial system-area

network, or a local-area network. A custom design is specific to a vendor

and is composed of networks of switches containing internal crossbars. An

example of commercial system-area network is Myrinet. Myrinet is a type of

local area network for high-performance interconnection of clusters \cite
{local area network for high-performance interconnection of clusters \cite}{%
BoCo95}. Ethernet is used as local-area network connection for networks of

workstation (NOW) architectures. In a network of workstations, each node of

the system is a workstation that collaborates with others via

message-passing.

In high-performance systems, efficient memory access schemes are needed.

Multiple levels caches are used to speed up data and instruction access.

Instruction reordering and data prefetching are used to avoid latency caused

by slow memories. This causes compilers to be left with the task of

generating efficient codes to take advantage of the hardware. Memory access

in shared memory systems has to incorporate cache coherence mechanisms to

avoid the access of non-valid data from cache.

Metacomputing is another approach to obtain high-performance systems. The

idea behind metacomputing is to obtain a highly powerful distributed system

composed of physically distributed computers to obtain the best resources

available to jointly solve a problem. The system should be transparent to

the users who will concentrate on the solution of their respective problems

and not on the computational requirements of the problem \cite{SmCa92}.

There are several issues in the development of such a global computing

platform, such as software compatibility, high-performance networks,

security, and user-friendly interfaces. Metacomputing involves the

interconnection of high-performance networks, implementing a distributed

file system, coordinating user access to different computational structures,

and making the environment easy-to-use and transparent to the user. There

are two projects addressing the metacomputing problem: Globus and Legion.

Globus is a project of Argonne National Laboratory and the University of

Southern California's Information Sciences Institute. The Globus\ project

intends to create the infrastructure to support the interaction of

computational systems physically located in geographically distributed

places. The goal of Globus is to understand application requirements for a

usable large-scale high-performance distributed computing environment \cite
{usable large-scale high-performance distributed computing environment \cite}%
{FoKe98}. There are three main research activities in Globus: (1) the

development of basic technology and high-level tools; (2) the construction

of a large-scale testbed; and (3) the execution of a real application to

evaluate the previous two activities.

The main component of the Globus system is the Globus Metacomputing Toolkit.

This toolkit is composed of the services needed for resource management,

communication, security, information, health and status, remote data access,

and executable management \cite{FoKe98}. Future services that will be

incorporated in the toolkit are accounting, auditing, and instrumentation.

The initial technologies for Globus have been implemented in the Globus

Ubiquitous Supercomputing Testbed (GUSTO). At the beginning of 1998, GUSTO

connected 17 sites in three different countries, 330 computers, and 3600

processors \cite{FoKe98}.

The Legion Project was developed at University of Virginia. Its goal is to

develop a worldwide virtual very powerful computer \cite{GrWu97}. The design

goals of Legion are scalability, usability, security, site autonomy, and

efficiency. Legion is an object-oriented system, that is, every component of

the system is represented by an object. This is the only programming

paradigm it supports. This is a main difference between Globus and Legion

since Globus supports different programming paradigms.

Legion should be responsible for protecting the user, the data, and the

resources from security problems. It should also manage data movement, fault

tolerance, and scheduling. Research areas involved in the project are:

parallel compilers, distributed systems software, object-oriented paradigms,

cryptography and cryptographic protocols. In the summer of 1995 the first

prototype of Legion was presented: the Campus Wide Virtual Computer (CWVC).

It consists of more than 100 workstations and a 18-processor IBM SP2

supercomputer in six buildings at the University of Virginia \cite{GrWu97}.

It posed some of the challenges of the nationwide Legion project to the

researchers.

\section{Environments and Tools}

We are proposing to improve integrated problem solving environments for the

development and efficient implementation of signal processing algorithms. We

will address the problem of raising the level of abstraction of the

information provided by tools to the signal processing specialist. In this

section we present some background on problem solving environments,

modeling, and high-performance computing tools.

\subsection{Problem Solving Environments\label{s:pse}}

A problem solving environment (PSE) is an integrated tool to provide

facilities for solving a class of problems in a convenient and efficient way

\cite{RiBo96}. It\ makes the computer system or architecture transparent to

the user, providing features which allow to solve problems ranging from

simple to complex. PSEs are designed with two main goals \cite{GaRi94}: they

should allow more people to solve more problems more quickly, and they

should enable users to solve otherwise unsolvable problems. Additionally,

PSEs should, in most cases, produce better and more reliable results.

Sekhar \emph{et al.} have developed a problem solving environment for signal

processing based on a dedicated artificial intelligence tool called OCAPI

\cite{ShMo94}. This environment brings the expertise of the signal

processing specialist to the signal processing user to solve an application.

This PSE is designed specifically for the solution of problems with the use

of current algorithms. This project does not address the design of signal

processing algorithms.

Another example of a problem solving environment is MATLAB$^{\registered }$.

MATLAB has become the dominant high-level language for mathematical and

scientific calculations. Six main reasons have been identified by Trefethen

\emph{et al.} for MATLAB's success: ease of use, high quality of numerical

algorithms, powerful graphic capabilities, high level of abstraction,

extensibility, and the use of matrices as basic data type [49-50].\nocite
{extensibility, and the use of matrices as basic data type [49-50].\nocite}{%
TrMe96} \nocite{MeTr97}

Users are able to use MATLAB for rapid development and testing of

algorithms. However, since MATLAB is an interpreted language, its

performance is not the same as the one obtained from compiled high level

languages such as C or Fortran. Several high performance versions of MATLAB

have been developed to improve its performance, two of which we will

describe: MultiMATLAB and RTExpress.

MultiMATLAB was developed at the Cornell Theory Center. It is an extension

of MATLAB to run in a distributed memory environment. In MultiMATLAB each

processor in the platform runs a MATLAB process which can communicate with

another via an interface module designed using MPI. Two different platforms

are supported: IBM SP2 and a network of workstations. MultiMATLAB's

performance is acceptable for coarse-grain distributed computing \cite
{performance is acceptable for coarse-grain distributed computing \cite}{%
MeTr97}, but for fine grain parallelism, its communication overhead is very

high when compared to directly programming in C and MPI. So there is a

trade-off between the rapid development of parallel algorithms and

performance \cite{Sa97}. Another feature of MultiMATLAB\ is that it supports

three different parallel programming paradigms: master/slave, single program

multiple data, and distributed matrices. This implies that the user has to

design an actual parallel algorithm to solve his/her problem. RTExpress has

another approach for high performance computing using MATLAB.

RTExpress$^{TM}$ is a compiler that converts MATLAB code into C with MPI for

high performance computers \cite{RTE}. It was developed by Integrated

Sensors, Inc. RTExpress supports the following platforms: Mercury, CPSI, IBM

SP2, Paragon, and a network of Sun workstations. RTExpress provides all

capabilities to compile, run, and execute MATLAB code on a target parallel

computer without having to recode the application. The following paradigms

are supported in\ RTExpress: data parallel, round robin, task parallel and

pipelining, and mixed mode. These paradigms are hidden to the user unless

he/she wants to specify compiler and/or run-time parameters.

\subsection{Performance Modeling}

In order to select alternatives in the development of code, the cost of each

alternative should be identified and measured. Performance models are used

for this purpose. Models are used to compare the efficiency of different

algorithms before doing an implementation. They can also be used to show

where improvements in an implementation can be made \cite{Fost95}.

A model should not be too simple so that estimates are far from reality or

too complex that it is impractical to use. An intermediate level of detail

is required to be useful to the programmer. It should be as simple but as

accurate as possible.

Performance models specify a metric as a function of algorithm or hardware

characteristic. For sequential architectures the metrics typically used are

execution time and memory utilization. Some of the metrics used for parallel

programs are execution time, speedup, and efficiency.

One of the most important metrics for parallel programs is execution time

\cite{KaFl90}. Execution time is the time elapsed between the execution of

the first instruction until the execution of the last instruction. Current

models for message passing environments describe execution time as a

combination of computation time, communication time, and idle time.

Efficiency measures the fraction of time processors spend doing useful work.

It is defined as

\end{document}
