VISPACKRoss T. WhitakerThe Scientific Computing and Imaging InstituteUniversity of UtahContact information:E-Mail: whitaker@cs.utah.eduWWW: http://www.cs.utah.edu/ whitakerUniversity of Utah, School of Computing, Technical Report: UUCS 08-0011Acknowledgments – This work is supported by the Office of Naval Research grants N00014-97-0227 andN00014-01-10033 and National Science Foundation grant CCR0092065. Input from many students, colleaguesand collaborators, including: Josh Cates and Samuel Burgiss Jr.AbstractVISPACK (volume-image-surface package) is a C++ library that includes matrix, image, and volumeobjects and tools for manipulating level-set surface models.2Contents1 Introduction 62 Design Attributes 62.1 Copy on Write or Data Handles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62.2 Direct Data Access . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82.3 Use of Templating . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92.4 Use of Standard File Formats . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102.5 Floating Point File Format . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102.6 Operator Overloading . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102.7 Image and Volume Processing Functions and Algorithms . . . . . . . . . . . . . . . . . . . . . 103 VISMatrix Library 103.1 Functionality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113.1.1 Basic Functionality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113.1.2 The VISVector Class . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123.1.3 Use of LAPACK . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123.2 Reference for the VISMatrix Library . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133.2.1 Interacting with the Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143.2.2 Overloaded VISMatrix/VISVector Operators . . . . . . . . . . . . . . . . . . . . . . . 173.2.3 General Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193.2.4 VISMatrix Specific Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 203.2.5 VISVector Specific Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 213.2.6 Creating VISVectors from Matrices and Vise-a-Versa . . . . . . . . . . . . . . . . . . . 214 Image Library 224.1 Functionality and Implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 224.1.1 Channels in Images . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2234.1.2 General Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 234.1.3 Image Processing Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 244.2 Reference for Image Library . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 244.2.1 Interacting with the Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 244.2.2 Overloaded Image Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 274.2.3 Math Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 294.2.4 General Image Processing Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 304.2.5 Image Processing Algorithms . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 325 Volume Library 335.1 Functionality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 335.1.1 General Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 345.1.2 Volume Processing Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 345.2 Reference for Volume Library . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 345.2.1 Interacting with the Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 355.2.2 Overloaded Volume Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 375.2.3 Math Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 395.2.4 General Volume Processing Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . 405.2.5 Volume Processing Algorithms . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 416 Util Library 416.1 Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 416.2 Reference for Util Library . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 426.2.1 Array Object . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 426.2.2 List Object . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 446.2.3 Math Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 447 Level-Set Surface-Modeling Library 4547.1 Surface Deformation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 467.2 Structure and Philosophy of the LSSM Library . . . . . . . . . . . . . . . . . . . . . . . . . . 477.3 The LevelSetModel Object . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 477.4 The VoxModel Object . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 487.5 Example: 3D Shape Metamorphasis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 487.6 Morph.h . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 497.7 Morph.C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 507.8 Main.C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52A Parameter File Objects 55A.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55A.2 Implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55A.3 Parameter file format . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55A.4 The ParameterFile object . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56B Code Examples 59B.1 matrix example . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61B.2 shotnoise example . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66B.3 fractal example . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6751 IntroductionVISPACK consists of five libraries: VISMatrix, Image, Volume, Util, and Level-Set Surface-Modeling. Theselibraries can be used separately or together when creating applications. Because VISPACK itself is dividedinto libraries we divide our discussion of VISPACK in this same way. Not only is this division natural forthe libraries (e.g., a matrix is much different from an image and they should not be included in the samelibrary), but this division also means that users need only understand the workings of the particular librarythey wish to use. We divide this document according to the libraries so that users can quickly find neededinformation.This format results in some redundancy, because there are similarities between the libraries. However,we have written this document so that these redundancies may easily be overlooked. For instance, accessmethods to ordered data structures (e.g. array or matrix) have similar conventions, but the descriptions ofthese methods are repeated for each object so that this document can be used as a reference.2 Design AttributesWe set forth three design goals for this library: usability, efficiency, and functionality. The philosophy behindthe usability goals is that the library should provide an environment in which users can easily implementsolutions to real problems. Our efficency goals are centered around the idea that the user should have accessto all the available resources of his/her machine. Our functionality goals ensure thoroughness in the library’sfunctionality. Table 1 shows our goals and the respective design attributes that support each goal.Table 1 indicates eight design attributes. In this section we give a discussions of these attributes. Theattributes presented are copy on write, direct data access, use of templating, use of standard file formats,floating point file format, operator overloading, and image and volume processing functions and algorithms.2.1 Copy on Write or Data HandlesWe use the idea of copy on write, or data handles, to make VISPACK more efficient for large data sets. Theadvantage of using this technique is that handles can be passed to and from functions rather than coping allTable 1: Goals and Design Attributes of VISPACKCategory Goal Design Attributeusability quick creation of test applications overloaded operators,copy on writecreation of large stable applications copy on write, direct data accessthat maintain effiencyeasy and accurate saving/loading of data use of standard file formatsfloating point image filesefficiency conservation of memory copy on write, use of templatesefficient implementation of multi-access functions copy on writelimit speed to that of machine direct data accessfunctionality basic matrix, image, and volume overloaded operatorsimage and volume processing implement functions and algorithmsability to create arrays ofmatricies, images, and volumes templating6the data for an object to pass to and from functions. Thus, copies (through either the copy constructor orassignment operator) are can be treated as deep copies which are deferred until they become necessary bythe copy-on-write protocol.To better understand this we examine an example. This example is the addition of three matrices and thenthe assignment of that result to a final matrix. Although this example uses matrices, we use this techniquefor images and volumes as well. An efficient way to implement the addition is the following 1:VISMatrix A(3,3),B(3,3),C(3,3);A=1.0f; //every element in A is 1B=2.0f; //every element in B is 2C=3.0f; //every element in C is 3VISMatrix ans;for (int r=0;rimage array;). An array of floats or ints is just as easily created. This templating makes the array andlist classes very powerful. Future releases of VISPACK will likely replace these objects with those fromthe Standard Template Library (STL).2.4 Use of Standard File FormatsWhen available we use standard file formats. We choose formats that are well known and have publiclyavailable libraries that can be distributed with our libraries. The matrix library uses a simple text format.The image library uses TIFF and FITS file formats. Because no standard format exists for saving volumesof data we do use a raw file format.2.5 Floating Point File FormatWe use a floating point file format in order to preserve information. For example, suppose a user wishesto apply a DOG (derivative of Gaussian) filter to an image and save the result. After, another applicationopens the saved result and performs further processing. The result of applying the DOG filter most likelycontains floating point values. Unless these floating point values are saved, the next processing step willbegin with inaccurate data. A floating point file format ensures that nothing is lost when saving files duringa multi-stage processing. In addition some information is inherently floating point, e.g., range data.2.6 Operator OverloadingUsed properly operator overloading gives users a convenient way to execute operations on an object. Usedimproperly, operator overloading can cause confusion and lead to poor readability. We, therefore, overloadonly basic operators, that is, operators with commonly known definitions. In addition, we clearly spell outthe functionality of all overloaded operators via inline documentation and in the Reference sections of thisdocument (sections 3.2, 4.2, 5.2, 7).2.7 Image and Volume Processing Functions and AlgorithmsThis library includes a number of general processing algorithms (e.g. edge detection) for images and volumes.These general processing functions provide a basis upon which algorithms may easily be built. In additionwe provide several specific algorithms for the image and volume objects. These specific algorithms allowusers to compare their results to the results of well-known algorithms we have implemented. These specificalgorithms are also useful as steps in a multi-stage process.3 VISMatrix LibraryThe ML (matrix library) contains a matrix class and a vector class. After reading the header file for theseclasses, anyone familiar with C++ should be able to quickly implement an application to test an idea thatrequires matrix manipulation. After more work with the library, users should begin to see how to efficiently10implement large scale, stable applications. Beyond this, an experienced programmer can, through directaccess to data buffers, write efficient implementations.The ML possesses many features common in linear algebra libraries and software packages. Because manyusers are already familiar with such systems, we do not discuss in detail the basic functionality of the ML.The purpose of this section is to briefly present the basic functionality of the ML and discuss the details ofthe implementation of the ML in a manner that will allow users to see how to get the most from the featuresof the ML.A few notes on convention:• In the ML we follow the convention of many linear algebra texts in that we refer to the size of andelements of a matrix in a row-by-column fashion (e.g. 4× 1 indicates four rows and one column).• Indexing in the ML begins at 0 (e.g. the first element of a vector is element number 0).Include files needed to use the library are:• matrix.h3.1 FunctionalityMany aspects of the functionality of the ML may be familiar to users. The following is a list of somecharacteristics of the ML:• Copy on write• Direct data access• Easy load/save functions• Use of LAPACKIn this section we first point users to sources of information concerning the basic functionality of the library.After, we discuss the functionality of the VISVector class. We present the VISVector class because theimplementation of this class my not be obvious to the user. Finally we describe functions in the library thatuse LAPACK (Linear Algebra PACKage). (We use LAPACK to implement several functions in the ML.)3.1.1 Basic FunctionalityThe two main sources for information on the basic functionality of the ML are the matrix.h file and section3.2 of this document. By examining these two sources of information, an experienced C++ programmershould be able to quickly learn about the VISMatrix and VISVector classes. The main function categoriesof basic functionality are listed in the table of contents in section 3.2.113.1.2 The VISVector ClassThe VISVector is a subclass of VISMatrix: vectors are n×1 matrices. To accomplish this we use inheritance.Because the VISMatrix class contains functions that do not apply to a vector we use private inheritance [7].An example of one such function is the svd function, which computes the singular value decomposition.Many computer vision applications work with sets of 2D or 3D points. Although the ML can be used inany way to work with sets of points, the ML is designed with a natural way of working with these points asvectors. In the ML, vectors are implemented as column matrices, specifically, a vector is actually a matrixwith 1 column. For instance, in the ML, multiplying a 3 element vector by a 3x3 matrix is illegal (run-timeerror). The ML sees this as multiplying a 3x1 matrix by a 3x3 matrix. The 3x3 matrix must be multipliedby the vector. The following code shows this:VISMatrix M(3,3);VISVector V(3);VISVector ans;ans = V*M; //invalidans = M*V; //correctIn addition, a natural way to store and recall sets of 2D or 3D points exists in the ML. VISVectors are storedin a matrix by creating a matrix with the number of columns equal to the number of points stored and thenumber of rows equal to the number of elements in the vectors. Two functions are used in this storing andrecalling process: VISMatrix concat(const VISVector& second) and VISVector vec(int col) const.(Other concat functions exist, see section 3.2.6.) The concat function stores the point (second) in thematrix (first), and the vec function recalls points from the matrix. The concat function can be used onan empty matrix; the result is simply a matrix that contains only the vector second. This is shown in theexample below:VISMatrix M;VISVector pt3D(3);for (int i=0;i<10;i++){pt3D=(float)i; //3D point is (i,i,i)M=M.concat(pt3D); //store the point}for (i=0;i<10;i++);{pt3D=M.vec(i); //recall the pointcout << "Here is point number " << i << endl << pt3D;}3.1.3 Use of LAPACKLAPACK (Linear Algebra PACKage) provides routines for solving systems of simultaneous linear equations,least-squares solutions of linear systems of equations, eigenvalue problems, and singular value problems [1].The associated matrix factorizations (LU, Cholesky, QR, SVD, Schur) are also provided. We use LAPACKroutines in the ML. This section describes how LAPACK is used in the ML and how to set up LAPACK inthe ML.The following functions of the ML (VISMatrix object) only are dependent on LAPACK (i.e. they will notwork without LAPACK):12• svd – singular value decomposition• det – determinant• QR – QR decomposition• inv, inverseSVD – inverse using SVD• inverseGJ – inverse using gauss jordan elimination• solvLinearEqn– solve a system of linear equationsWe create the above functions using functions that convert a VISMatrix to/from LAPACK’s data formats.This allows the use of any LAPACK function. For example the svd function. We convert the matrix beingprocessed to the LAPACK format. We then apply the LAPACK function sgesvd . Finally we convert theresulting LAPACK matrix back to a VISMatrix object and return that object as the result.The conversion functions metioned above are available for users. For example, suppose a user wishes toimplement an algorithm that requires determining the eigenvalues and eigenvectors of a matrix, which doesnot currently exist in the ML. A user first inherits the matrix class and then adds a function to the inheritedclass that utilizes the sgeev () function in the C library of LAPACK to find eigenvalues and eigenvectors.The function of the inherited class would:1. Convert the inherited matrix into a LAPACK format matrix2. Apply the sgeev () LAPACK function to the LAPACK format matrix3. Convert the result back to the inherited matrix formatThe inheritance is necessary because the conversion functions are protected members of the ML. The func-tions are protected rather than public because users should add basic functionality as memeber functionsand then use this basic functionallity to write functions outside the class.The functions available for conversion to/from LAPACK format are: float* VISMatrix::createLAPK()const,void VISMatrix::becomeLAPK(float* lapk, int r, int c),float* VISVector::createLAPK() const,and void VISVector::becomeLAPK(float* lapk, int n). The createLAPK functions create a LAPACKmatrix/vector from an ML matrix/vector respectively. The becomeLAPK functions assign the values from aLAPACK matrix/vector to a ML matrix/vector respectivelly. See section 3.2 for more information on thesefunctions.We include LAPACK with this release of the library. Users should downlad LAPACK updates and docu-mentation on the World-Wide Web at http://www.netlib.org/lapack/ .3.2 Reference for the VISMatrix LibraryThis reference section lists the functions of the ML. The functions are listed according to category. If onlyone function declaration is listed the function has the same declaration in the both the VISMatrix andVISVector classes. Otherwise, all the declarations are listed.13Table of contents:3.2.1 Interfacing with the dataPeeking and PokingDirect Access to the DataConstructing a VISMatrix/VISVectorFile I/OConverting to/from LAPACK Format3.2.2 Overloaded VISMatrix/VISVector operators3.2.3 General Functions3.2.4 VISMatrix Specific Functions3.2.5 VISVector Specific Functions3.2.6 Creating VISVectors from VISMatrixices and Vise-a-Versa3.2.1 Interacting with the DataThis section includes several subsections entitled: Peeking and Poking, Direct Access to the Data, Con-structing a matrix, File I/O, and Converting to/from LAPACK Format. Peeking and Poking lists the mainfunctions for accessing data and changing data of a matrix/vector. Direct access to the Data Buffers listsfunctions that return pointers to the memory block holding the data of a matrix. Constructing a matrixlists ways to construct a matrix, and File I/O lists methods for inputing/outputing matricies and vectorsfrom/to files.Peeking and Pokingpoke – float& VISMatrix::poke(int row,int col)float& VISVector::poke(int n)Assign a value to the element of a matrix/vector indicated by row and col (or n). This function isalways used on the lhs of an assignment.Example:VISMatrix M(4,4);M.poke(2,3)=1.44;peek -- float VISMatrix::peek(int row,int col)float VISVector::peek(int n)Returns the value of the element of a matrix/vector indicated by row and col (or n). Example:VISMatrix M(4,4);M.poke(2,3)=1.44;float m=M.peek(2,3);operator() – float VISMatrix::operator()(int row,int col)float VISVector::operator()(int n)Same as peek but easier to type. Returns the value of the element of a matrix/vector indicated by rowand col (or n). Example:VISMatrix M(4,4);M.poke(2,3)=1.44;float m=M(2,3);VISVector V(2);V.poke(2)=1.55;float v=V(2);14peekROI – VISMatrix peekROI(int startrow, int endrow, int startcol, int endcol) constVISVector peekROI(int startelement,int endelement) constThese functions returns a region of interest from a matrix/vector. The region is defined by startrow,endrow and startcol, endcol OR startelement, endelement.putROI – pokeROI(int row,int col,VISMatrix roi)pokeROI(int row,int col,VISVector roi)pokeROI(int row,VISVector roi)These function allows a user to place a region of interest into a matrix/vector. The region itself iseither a matrix or a vector. The upper left corner of the region of interest is indicated by row and col.Direct Access to the Data BuffersdataPtr – const float* VISMatrix::dataPtr(int i = 0) constconst float* VISVector::dataPtr(int i = 0) constThis function returns a pointer to a block of float data representing a matrix/vector. The functionreturns, by default, the pointer to the begining of the data block, otherwise, the pointer at the positionindicated by i is returned. The matrix data is stored in row major order.data float data(int i = 0) constThis function returns the floating point value at location i in the data block if i is specified. Otherwise,the value at the first position in the data block is returned. We store the matrix data in row majororder.Constructing a VISMatrix These functions are listed according to the input of the function.default – VISMatrix::VISMatrix()VISVector::VISVector()An empty matrix/vector is created.other – VISMatrix::VISMatrix(VISMatrix& other)VISVector::VISVector(VISVector& other)These functions create a matrix/vector from another matrix/vector. The matrix/vector created is thesame other.row,col/n – VISMatrix::VISMatrix(int row, int col)VISVector::VISVector(int n)These functions create a matrix/vector of a specified size. The size of a matrix is rowxcol. The sizeof a vector is n.char* filename – VISMatrix::VISMatrix(char* filename)VISVector::VISVector(char* filename)These functions create a matrix/vector from a file. The file must contain only one matrix/vector wherethe output format for a VISMatrix is:numberofrows numberofcolsM0,0 M0,1 M0,2 . . .M1,0 M1,1 M1,2 . . .. . . . . .and the output format for a VISVector is:numberofelements15V0V1V2. . .Example:VISMatrix M("file.mat");File I/Ooperator<< – ostream& operator<<(ostream& os, const VISMatrix& m)ostream& operator<<(ostream& os, const VISVector& m)Print a matrix to any ostream. Output format for a VISMatrix is:numberofrows numberofcolsM0,0 M0,1 M0,2 . . .M1,0 M1,1 M1,2 . . .. . . . . .Output format for VISVector is:numberofelementsV0V1V2. . .Example:VISMatrix M(3,2);M=2.0f;cout << M;ofstream file("output");file << M;operator>> – istream& operator>>(istream &is, VISMatrix& m)istream& operator>>(istream &is, VISVector& m)Read a matrix to any ifstream. Input format is the same as output format (see operator<<). Example:VISMatrix M(3,2);M=2.0f;ifstream file("input");file >> M;Converting to/from LAPACK FormatbecomeLAPK – void VISMatrix::becomeLAPK(float* lapk, int r, int c),void VISVector::becomeLAPK(float* lapk, int n).These functions assign the values from a LAPACK format matrix/vector to an ML VISMatrix/VISVectorrespectively.16createLAPK – float* VISMatrix::createLAPK() const,float* VISVector::createLAPK() const,and These functions create a LAPACK format matrix/vector from an ML VISMatrix/VISVectorrespectively.3.2.2 Overloaded VISMatrix/VISVector OperatorsSome of these functions may have the same name but different inputs and outputs. In these cases the labelother means an input argument of type VISMatrix or VISVector; whereas, the label s refers to a floatingpoint number.operator=(other) – VISMatrix& VISMatrix::operator=(const VISMatrix& other)VISVector& VISVector::operator=(const VISVector& other)VISVector& VISVector::operator=(const VISMatrix& other)Assigns one matrix/vector to another. This function does not copy the data. (See section 2.1.) Amatrix may be assigned to a vector only if the matrix has 1 column. Example:VISMatrix M1(4,4);M1=1.0f;VISMatrix M2;M2=M1;VISMatrix M3(4,1);M3=5.0f;VISVector V1;V1=M1; //invalidV1=M3; //o.k.operator=(s) – VISMatrix& VISMatrix::operator=(float s)VISVector& VISVector::operator=(float s)Assigns each element in a matrix/vector to the value of the float s. Example:VISMatrix A(4,4);A=3.141592;//Every element of A is equal to 3.141592operator+(other) – VISMatrix VISMatrix::operator+(const VISMatrix& other) constVISVector VISVector::operator+(const VISVector& other)Adds two matrices/vectors. The matrices/vectors must be the same size. Example:VISMatrix A(4,4),B(4,4),C;A=1.0f;B=2.0f;C=A+B;//now C=3.0foperator+(s) – VISMatrix VISMatrix::operator+(float s) constVISVector VISVector::operator+(float s) constAdds the value of s to a matrix/vector.Example:VISMatrix A(4,4),B;A=1.0f;B=A+2.0f;//now B=3.0f;17operator-(other) – VISMatrix VISMatrix::operator-(const VISMatrix& other) constVISVector VISVector::operator-(const VISVector& other) constSubtracts two matrices/vectors. The matrices/vectors must be the same size.Example:VISMatrix A(4,4),B(4,4),C;A=3.0f;B=1.0f;C=A-B;//now C=2.0foperator-(s) – VISMatrix VISMatrix::operator-(float s) constVISVector VISVector::operator-(float s) constSubtracts the value of s from a matrix/vector. Example:VISMatrix A(4,4),B;A=10.0f;B=A-2.0f;//now B=8.0f;operator-() – VISMatrix VISMatrix::operator-() constVISVector VISVector::operator-() constNegates a matrix/vector. Example:VISMatrix A(4,4),B;A=3.0f;B=-A;//now B=-3.0foperator*(other) – VISMatrix VISMatrix::operator*(const VISMatrix& other) constVISMatrix VISVector::operator*(const VISMatrix& other) constVISVector VISMatrix::operator*(const VISVector& other) constMultiplies a matrix/vector by a matrix/vector (other). Examine the example of A=B*C. The number ofcolumns of the B must be equal to the number of rows of C (other). Implementing vectors as matriceswith only one column, as we have done, means that if B is a vector and C a matrix the result, A, is amatrix. On the other hand if B is a matrix and C is a vector the result, A, is a vector. Dot product andcross product operations are provided in the vector class. Example:VISMatrix A(4,4),B(4,4),C;A=3.0f;B=4.0f;C=A*B;//now C=48.0f and is 4x4C=B*A;//still C=48.0f, because A and B are square A*B=B*AVISVector D(4),E;D=2.0f;E=D*B;//error!E=B*D;//now E=24.0f and has four elements (4x1)VISVector FF=E*A;//error!F=A*E;//now F=288.0f and has four elements (4x1)operator*(s) – VISMatrix VISMatrix::operator*(float s) constVISMatrix operator*(float s,const VISMatrix m)VISVector VISVector::operator*(float s) constVISVector operator*(float s,const VISVector v)Multiplies each value in a matrix/vector by s. The non-member functions allow the order of multipli-cation be reversed (i.e. M ∗ 2.0f or 2.0f ∗M). Example:VISMatrix A(4,4),B;A=2.0f;B=A*2.0f;//now B=4.0f;18B=4.0f*A;//now B=8.0f;operator/(s) – VISMatrix VISMatrix::operator/(float s) constVISMatrix operator/(float s,const VISMatrix m)VISVector VISVector::operator/(float s) constVISVector operator/(float s,const VISVector v)Divides each value in a matrix/vector by s or s by each value in a matrix/vector.Example:VISMatrix A(4,4),B;A=2.0f;B=A/4.0f;//now B=0.5f;B=8.0f/A;//now B=4.0f;operator==(other) – VISMatrix VISMatrix::operator==(const VISMatrix& other) constVISVector VISVector::operator==(const VISVector& other) constWhen used as A==B this function returns a matrix/vector containing 0/1’s indicating equality/inequalityrespectively of each element of A to the respective element of B. Example:VISMatrix A(4,4),B(4,4),C;A=2.0f;B=2.0f;A.poke(0,0)=1.0f;C = A==B;//now all elements of C=1.0f except 0,0 which is 0.0f;To determine if every element of A is equal to the respective elements in B use the following:if((A==B).min()) ... //tests if all elements of A are equal to Boperator==(s) – VISMatrix VISMatrix::operator==(float s) constVISVector VISVector::operator==(float s) constReturns a matrix/vector of 0/1’s indicating equality/inequality respectively of each element to thevalue of s. Example:VISMatrix A(4,4),B;A=2.0f;A.poke(0,0)=1.0f;B = A==2.0f;//now all elements of B=1.0f except 0,0 which is 0.0f;To determine if every element of A is equal to s use the following:if ((A==2.0f).min()) //tests if all elements of A are equal to k... //where k is a floating point value3.2.3 General Functionst – VISMatrix VISMatrix::t() constVISMatrix VISVector::t()This operator returns the transpose of a matrix/vector. Because a vector has only 1 column, thetranspose of a VISVector is a VISMatrix that has 1 row. Example:VISVector A(4),B;A=2.0f;//A has 1 column (4x1)B=A.t();//B is the transpose of A, B has 1 row (1x4)max – float VISMatrix::max() constfloat VISVector::max() constfloat VISMatrix::max(int& r,int& c) const19float VISVector::max(int& n) constThese functions return the maximum value in a matrix/vector. In addition, in the case of a matrix,this function assigns the row and column of the location of the maximum value to r and c. In the caseof a vector the function assigns the row of the location of the maximum value to n.min – float VISMatrix::min() constfloat VISVector::min() constfloat VISMatrix::min(int& r,int& c)float VISVector::min(int& n)These functions return the minimum value in a matrix/vector. In addition, in the case of a matrix,this function assigns the row and column of the location of the minimum value to r and c. In the caseof a vector the function assigns the row of the location of the minimum value to n.mean – VISVector VISMatrix::mean() constVISVector VISMatrix::meanOfCols() constfloat VISVector::mean() constThe function VISMatrix::mean returns the mean values of each row in a matrix. The functionmeanOfCols returns the mean values of each column in a matrix. The function VISVector::meanreturns the mean of a vector.3.2.4 VISMatrix Specific Functionstranspose – See function t in section 3.2.3.inv – VISMatrix VISMatrix::inv() constVISMatrix inverseSVD() constVISMatrix inverseGJ() constThe first two functions compute the inverse of a matrix using SVD (singular value decomposition).The second function computes the inverse of a matrix using Gauss Jordan elimination. These functionswill not work without LAPACK. Example:VISMatrix A(4,4),B;A=2.0f;B=A.inv();//B is the inverse of Anorm – float VISMatrix::normOfCol(int col = 0) constfloat VISMatrix::normOfRow(int row = 0) constThe function normOfCol returns the 2-norm of a column of a matrix. The function normOfRow returnsthe 2-norm of a row of a matrix. Default is column or row number is 0.det – float VISMatrix::det() constThis function finds the determinant of a matrix using SVD. This function will not work withoutLAPACK.svd – void svd(VISMatrix &U, VISMatrix &W, VISMatrix &V) constThis function performs singular value decomposition. The matrices from the decomposition are placedin the matrices U, W, and V that are inputs to the function. This function will not work withoutLAPACK.i – VISMatrix i(int size)Returns an identity matrix that is has size rows and size columns.cov – VISMatrix VISMatrix::cov() constFor matrices where each row is an observation and each column a variable this function returns thecovariance matrix.20QR – void QR(VISMatrix& Q, VISMatrix& R) constReturns the QR decomposition of a matrix.solvLinearEqn – VISVector solvLinearEqn(const VISMatrix& A,const VISVector& b)Solves the linear equations:Ax = b, (1)where A is a square matrix, b and x are vectors. The returned result is represented by x in thisequation.3.2.5 VISVector Specific Functionsnorm – float VISVector::norm() constThe function norm returns the 2-norm of a vector.cross – VISVector VISVector::cross(const VISVector& other) constThis function computes the cross product of two vectors.dot(other) – float dot(const VISVector& other) constThis function returns the dot product of two vectors.3.2.6 Creating VISVectors from Matrices and Vise-a-Versavec – VISVector VISMatrix::vec(int col = 0) constVISVector VISMatrix::vecfromcol(int col = 0) constVISVector VISMatrix::vecFromRow(int row = 0) constVISVector VISMatrix::vecFromDiag() constThese functions create a vector from a matrix. The function vec and vecfromcol create a vectorfrom a column of a matrix. (These functions are the same.) The default column is 0. The functionvecFromRow creates a vector from a row of a matrix. The default row is 0. The function vecFromDiagcreates a vector from the diagonal of a matrix.concat – VISMatrix concat(const VISMatrix& second)VISMatrix concat(const VISVector& second)These functions concatenate the columns of a matrix/vector to the columns of a second matrix/vector.Recall that a vector is treated as a matrix having one column. These functions will work even if oneof the matrices/vectors is empty; the matrix/vector containing data is returned by the function.Example:Suppose A =[1 2 34 5 6]and B =[7 8 910 11 12]then A.concat(B) =[1 2 3 7 8 94 5 6 10 11 12].concatRow – VISMatrix concatRow(const VISMatrix& second)VISVector concatRow(const VISVector& second)These functions concatenate the rows of a matrix/vector to the rows of a second matrix/vector. Be-cause a vector has only one column there are no functions to concatenate the rows of a vector tothe rows of a matrix or vise-a-versa. To concatenate a vector (as a row) onto the rows of a matrixuse M.concatRow(V.t()), where M is a matrix having n columns and V is a vector with n elements.Example:21Suppose A =[1 2 34 5 6]and B =[7 8 910 11 12]then A.concat(A) =1 2 37 8 94 5 610 11 12.4 Image LibraryThe IL (image library) contains two main classes: the VISImage class and the VISImageFile class. TheVISImageFile class is used only for file I/O and the VISImage class is used for all other functionality. Theheader files, image.h and imagefile.h, which include these two classes, provide much information. Theseheader files in conjunction with this document give users complete information.Our presentation of the IL includes a discussion of some of the functionality and implementation issues ofthe library and a reference section. Our discussion of functionality and implementation issues describes thegeneral functionality of the library and image processing functionality. The reference provided for the ILlists each function according to category.In this document we refer to the size of and pixels of an image (matrices are different) in a column-by-rowfashion (e.g. 1x4, four columns one row). In addition, the IL follows this same convention.Include files needed to use the library are:• image.h• imagefile.h• imageRGBA.h4.1 Functionality and ImplementationWe divide the functions of the IL into two categories: general functions and image processing functions.General functions include functions that allow the user to interact with data stored in an image, file I/Ofunctions, overloaded operators, and basic math functions. Image processing functions refer to functionsthat are more specific to image processing. Filtering, resampling, and distance transform functions are allconsidered image processing functions. In addition IL includes more complex image processing algorithms;for example, canny edge detection. We divide our discussion of the functionality of this library into the twocategories: general functions and image processing functions.4.1.1 Channels in ImagesImages can contain multiple 2D arrays of data. Each 2D array is called a channel. By default, all imagesare of a single channel and all access methods pertain to the first channel. However, channel numbers canbe given for constructors and access methods, and virtually all of the processing routines will apply to all ofthe channels. The use of channels is particularly useful when dealing with unpacked color images.224.1.2 General FunctionsIn this section we point users to sources of information that explain the specifics of the general functionsin the IL, and we discuss some attributes of the general functions of the IL that may not be obvious fromthe information given. Here, our discussion of general functions is broken into parts. We describe functionsthat allow access to the data in an image, file I/O functions, overloaded operators, math functions, imageprocessing functions, and miscellaneous functions. Beyond this information, detailed information for allthe functions is available from other sources. Descriptions of the general functions (excluding a file I/Ofunctions) can be found in the image.h header file. Detailed information about file I/O can be found in theimagefile.h header file. In addition section 4.2 of this document provides descriptions of each function inthe IL.Accessing data The IL uses copy on write as described in section 2.1. This means that there are twodifferent functions for accessing the data of a pixel in an image. One function for use on the left hand sideof an assignment (read/write) and another for use on the right hand side of an assignment (read only). Anassignment in IL looks like the following:VISImage A(128,128);VISImage B(128,128); A=1.0f; //every pixel in A is 1B=2.0f; //every pixel in B is 2A.at(4,4) = B(4,4); //now pixel 4,4 of A = pixel 4,4 of B, which = 1//the lhs uses at(column,row) and the rhs uses//operator(column,row)File I/O The IL uses a second class, VISImageFile to accomplish file I/O. Information about this classis contained in the file imagefile.h. Here we present a piece of example code that reads a floating pointimage from a file, adds a value to the image, and then saves the result in a floating point file format:VISIm im; // a typeless base classVISImage input; //input imageVISImage output; //output imageVISImageFile im file;if (!((im=im file.read("filename")).isValid())){fprintf(stderr,"Error reading image file %s.\ n","filename");exit(0);}input=VISImage(im);//Cast input image as a floating point imageoutput = input + 2.0f; //Add 2 to the input image assign result to outputfile.write fits(output,"out.fit");To read/write images of type other than float use functions read tiff and write tiff (see section 4.2.1).For more usage information see also the file imagetest/imagetest.cxx in the source code.Overloaded Operators The VISImage class contains overloaded operators. A value may be added to orsubtracted from an image using + and - respectively. An image may be multiplied by or divided by a valueby * and /, respectively. Applying +, -, *, and / to two images adds, subtracts, multiplies, and divides therespective pixels in the two images.23Math Functions Many math functions exist in the IL–for example: square root, natural logarithm,exponential, absolute value, etc. All of these functions have short function names that are designed to beeasily recalled.4.1.3 Image Processing FunctionsGeneral image processing functions exist that:• apply masks to an image• smooth an image• determine derivatives of an image• resample an image• perform distance transforms on an image• corupt an image with noiseDetails about these functions can be found in the image.h header file and in section 4.2.4 of this document.In addition to general image processing functions, several image processing algorithms exist in the IL.Information about image processing algorithms can also be found in the image.h header file or section 4.2.5.4.2 Reference for Image LibraryThis reference section lists the functions of the IL. The functions are listed according to category.Table of contents:4.2.1 Interacting with the dataPeeking and PokingDirect Access to the Data BuffersConstructing an ImageFile I/O4.2.2 Overloaded Image Operators4.2.3 Math Functions4.2.4 General Image Processing Functions4.2.5 Image Processing Algorithms4.2.1 Interacting with the DataThis section includes several subsections entitled: Peeking and Poking, Direct Access to the Data Buffers,Constructing an Image, and File I/O. All the functions listed here are members of the VISImage classunless otherwise noted. Peeking and Poking lists the main functions for accessing data and changing data ofan image. We describe the idea of direct data access in section 2.2; here the Direct Data Access descriptionlists functions related to this idea. Constructing an Image lists ways to construct an image, and File I/Olists methods for inputing/outputting images from/to files.24Peeking and Pokingpoke – T& poke(unsigned int x, unsigned int y) constT& poke(unsigned int x, unsigned int y, unsigned int ch) constAssign a value to the pixel of an image indicated by x (column), y (row), and ch (channel). Thisfunction is always used on the lhs of an assignment. Example:VISImage I(4,4);I.poke(2,3)=1.44;peek – const T& peek(unsigned int x, unsigned int y) constconst T& peek(unsigned int x, unsigned int y, unsigned int ch) constReturns a value of the pixel of an image indicated by x (column), y (row), and ch (channel). Example:VISImage I(4,4);I.poke(2,3)=1.44;float m=I.peek(2,3);operator() – const T& operator()(unsigned int x, unsigned int y) constconst T& operator()(unsigned int x, unsigned int y, unsigned int ch) constSame as peek. Returns a value of the pixel of an image indicated by x (column), y (row), and ch(channel). Example:VISImage I(4,4);I.poke(2,3)=1.44;float m=I(2,3);putChannel – void putChannel(VISImage& other)void putChannel(VISImage& other, unsigned int ch)These functions place a 2D image into a channel of an image. If ch is not specified (first function),a new channel is created, otherwise (second function), the 2D image is placed into channel ch of amulti-channel image.channel – VISImage channel(unsigned int ch)This function returns the 2D image at channel ch of a multi-channel image.interp – T interp(float x, float y, unsigned int ch) constT interp(float x, float y) constThese functions are similar to peek and operator(), but they interpolate between pixel locations ofthe image to determine the value at the location represented by x, y, and ch. These functions uselinear interpolation.interpNoBounds – T interpNoBounds(float x, float y, unsigned int ch) constT interpNoBounds(float x, float y) constThese two functions are the same as the interp functions except they do not do any bounds checking.By bounds checking we mean ensuring that the values of x, y, and ch are valid (i.e. within limits).Hence, these functions are faster than the previous two, but are not as safe to use.getROI – VISImage getROI(unsigned int x pos, unsigned int y pos,unsigned int w roi, unsigned int h roi)This function returns a region of interest from an image. The position (x pos,y pos) indicates theupper left corner of the image. In addition the region has size indicated by w roi,h roi– its width andheight.putROI – void putROI(const VISImage& image in, unsigned int x pos, unsigned int y pos)This function allows a user to place a region of interest into an image. The region itself is an imageand has a position (x pos, y pos)–its upper left corner.25printData – void printData() constThis does a floating point printf on all of the data. Be careful if you have a large image.Direct Access to the Data Buffersrep – const VISImageRep* rep(unsigned int ch) constconst VISImageRep* rep() constThis routine returns a const pointer to the data of a chosen channel of an image if ch is indicated,otherwise, the function returns a const pointer to the data of the first channel of the image.repRef – VISImageRep* repRef(unsigned int ch)VISImageRep* repRef()This routine returns a non-const pointer to the data of a chosen channel of an image if ch is indicated,otherwise, the function returns a non-const pointer to the data of the first channel of the image. Thisroutine copies the data. (See section 2.1.)Constructing an Image Most of these functions are listed according to the input of the function.default – VISImage()An empty image is created.image – VISImage(const VISImage& image)This function creates an image from another image. The image created is exactly the same as image.w,h,ch – VISImage(unsigned int w, unsigned int h)VISImage(unsigned int w, unsigned int h, unsigned int ch)These functions create an image of a specified size, where w is width, h is height, and ch is the numberof channels.buffer VISImage(unsigned int w, unsigned int h, unsigned int ch, T** buf) This function cre-ates an image w x h with ch channels, using buf as the data. The pointer buf is a pointer to anarray of pointers to individual images (channels). Each image buffer is a w × h long 1D array of typeT.createToSize VISImage createToSize() constReturns an image of the same size of the image being operated on. Example:VISImage I1(128,128);VISImage I2;I2 = I1.createToSize; //I2 is now 128x128, values in the image are unknownFile I/O There is only one read function in the IL. This function automatically determines the file typeto read. File formats supported by the read function are discussed below. There are several write functionsfor the different file formats supported by the IL.read VISIm read(const char* fname)This function reads any file type supported by the IL. These file types are tiff and fits. The tiff fileformat supports most image types, including color, but does not support floating point. The fits fileformat supports all image types, but does not have a specific type for color. Color images in fits mustbe handled as multi-channel images. Example (reading a floating point fits format file):26VISImage input; //input imageVISImageFile file;if (!((input=file.read("filename.fit")).isValid())){fprintf(stderr,"Error reading image file %s.\ n","filename");exit(0);}input=VISImage(input);//Cast input image as a floating point imagewrite tiff – int write tiff(const VISImageRGBA& image, const char* fname)int write tiff(const VISImage& image, const char* fname)int write tiff(const VISImage& image, const char* fname)int write tiff(const VISImage& image, const char* fname)int write tiff(const VISImage& image, const char* fname)These functions write any image type as a tiff file. The file name is indicated by fname.write fits – int write fits(const VISImage& image, const char* fname)int write fits(const VISImage& image, const char* fname)int write fits(const VISImage& image, const char* fname)int write fits(const VISImage&, const char*)These functions write any image type, excluding color images, as a fits format file. The file is namedfname.write iv – int write iv(const VISImage& im, const char* filename)int write iv(const VISImage& x, const VISImage& y, const VISImage&z,const char* filename)int VISImageFile::write iv(const VISImage& im, const VISImageRGBA& im color, constchar* filename)These functions create a file that is a surface mesh in open inventor format. The mesh is tesslatedon the image grid (2 triangles per group of four pixels). In the first case the image coordinates aretreated as x and y positions and the value of the image at each coordinate is treated as z—this givesthe graph of an image. In the second case, the grid coordinates define the topology of the mesh, butthe x, y, and z positions for each image coordinate are given by the values of the 3 input images. Thethird case is the same as the first but also includes a color (material property) at each vertex from thevalues given by the color image. For multiple input images, they must be the same size.4.2.2 Overloaded Image OperatorsThe input to the functions listed here is either an image or a value. We list functions according to operationname and input type. The label image indicates an object of type VISImage, and the label value refersto a value. (The function name indicates the type of value.)operator=(image) – VISImage& operator=(const VISImage& from)Assigns one image to another. This function does not copy the data. (See section 2.1.) Example:VISImage A(128,128);A=1.0f; //every pixel in A is 1VISImage B;B=A; //B is the same as Aoperator=(value) – VISImage& operator=(T value)Sets each pixel in an image to value. Example:VISImage A(128,128);A=3.141592;//Every pixel of A is equal to 3.14159227operator+(image) – VISImage operator+(const VISImage& image) constAdds two images. The images must be the same size. Example:VISImage A(128,128),B(128,128),C;A=1.0f;B=2.0f;C=A+B;//now C=3.0foperator+(value) – VISImage operator+(const VISImage& from, T value)VISImage operator+(T value,const VISImage& from)Adds value to an image.Example:VISImage A(128,128),B;A=1.0f;B=A+2.0f;//now B=3.0fB=3.0f+A;//now B=4.0foperator-(image) – VISImage operator-(const VISImage& image) constSubtracts one image from another image. The images must be the same size. Example:VISImage A(128,128),B(128,128),C;A=1.0f;B=2.0f;C=B-A;//now C=1.0foperator-(value) – VISImage operator-(const VISImage&from, T value)VISImage operator-(T value,const VISImage& from)Subtracts a value from an image or subtracts an image from a value. Example:VISImage A(128,128),B;A=1.0f;B=2.0f-A;//now B=1.0fB=A-3.0f;//now B=-2.0foperator*(image) – VISImage operator*(const VISImage& image) constMultiplies each pixel in an image by the respective pixel in another image, image. Example:VISImage A(128,128),B(128,128),C;A=5.0f;B=2.0f;C=B*A;//now C=10.0foperator*(value) – VISImage operator*(const VISImage& from, T value)VISImage operator*(T value, const VISImage& from)Multiplies each pixel in an image by value. Example:VISImage A(128,128),B;A=2.0f;B=2.0f * A;//now B=4.0fB=A * -3.0f;//now B=-6.0foperator/(image) – VISImage operator/(const VISImage& image) constDivides each pixel in an image by the respective pixel in image. Example:VISImage A(128,128),B(128,128),C;A=5.0f;B=2.0f;C=B/A;//now C=0.4operator/(value) – VISImage operator/(const VISImage& from, T value)VISImage operator/(T value, const VISImage& from)28Divides each pixel in an image by value. Example:VISImage A(128,128),B;A=2.0f;B=10.0f/A;//now B=5.0fB=A / 4.0f;//now B=0.5operator+=(image) – VISImage& operator+=(const VISImage& image)Adds each pixel in an image to the respective pixel in image and assigns the result to the image beingoperated on.operator+=(value) – VISImage& operator+=(const T& value)Adds value to each pixel and assigns the result to the image being operated on.operator-=(image) – VISImage& operator-=(const VISImage& image)Subtracts each pixel in image from the respective pixel in the image being operated on and assigns theresult to the image being operated on.operator-=(value) – VISImage& operator-=(const T& value)Subtracts value from each pixel and assigns the result to the image being operated on.operator*=(image) – VISImage& operator*=(const VISImage& image)Multiplies each pixel by the respective pixel in image and assigns the result to the image being operatedon.operator*=(value) – VISImage& operator*=(const T& value)Multiplies each pixel by value and assigns the result to the image being operated on.operator/=(image) – VISImage& operator/=(const VISImage& image)Divides each pixel by the respective pixel in image and assigns the result to the image being operatedon.operator/=(value) – VISImage& operator/=(const T& value) Divides each pixel by value and as-signs the result to the image being operated on.4.2.3 Math Functionspower(int exponent) – VISImage power(const int& exponent) constThis function raises the image to a power.sqrt() – VISImage sqrt() constPerforms the square root of an image.div by(value,zeroCondition) – VISImage div by(T value, T zeroCondition) constDivides an image by value; if division by zero is encountered the pixel value in the returned image isequal to zeroCondition.div(image,zeroCondition) – VISImage div(const VISImage& image, T zeroCondition) constDivides the pixels in an image by the respective pixels in image; if division by zero is encountered thepixel value in the returned image is equal to zeroCondition.ln() – VISImage ln() constReturns natural log of the image.exp() – VISImage exp() constReturns the exponential of the image.29pos() – VISImage pos() constReplaces negative numbers with zeros. Does not change positive numbers.neg() – VISImage neg() constReplaces positive numbers with zeros. Does not change negative numbers.abs() – VISImage abs() constReturns the absolute value of the image.sign() – VISImage sign() constReplaces all zeros with +1. Replaces all positive numbers with +1. Replaces all negative numberswith -1.max() – T max() constReturns max value in image.min() – T min() constReturns min value in image.max(image) – VISImage max(const VISImage &other) constReturns an image containing the max of other and the image being operated on.min(image) – VISImage min(const VISImage &other) constReturns an image containing the min of other and the image being operated on.sum() – float sum() constReturns the sum of all pixels in the image.average() – float average()Returns the average of all pixels in the imagescale(value) – VISImage scale(float value) const;Multiplies each pixel by value.4.2.4 General Image Processing Functionsmedian – VISImage median(int window size) constMedian filter where window size indicates the size of the mask used.gaus – VISImage gauss(float sigma) constGaussian smoothing. The standard deviation of the Gaussian filter is sigma.convolve – VISImage convolve(const VISImage& kernel) constConvolution with any filter, which is given in the form of an image. The image library containsglobal functions for constructing images that correspond to Gaussian smoothing filters and discretederivatives.mask – VISImage mask(const VISImage& kernel) constDetermines the inner product of a mask and an area in an image.maskFloat – float maskFloat(const VISImage& mask, unsigned int center x,unsigned int center y, unsigned int x pos, unsigned int y pos) constArray* maskFloat(const Array< VISImage >& masks,unsigned int x pos, unsigned int y pos) constDetermines the inner product of a mask and an area in an image. In the first function, the center ofthe mask is indicated by (center x, center y) and the center of the mask is positioned at (x pos,30y pos) of the image and then the inner product is calculated. The second function does the same asthe first except an array of masks is applied to the image, and thus, an array of floats is returned.scaleToRGB – VISImage scaleToRGB()Maps the grey-scale values in an image to the interval [0,255] using a linear map.dx – VISImage dx() constVISImage dx(unsigned int order) constThese functions determine the derivative in x direction (along rows). The second function takes theorder argument indicating the order of the derivative to be taken.dy – VISImage dy() constVISImage dy(unsigned int order) constThese functions determine the derivative in y direction (along columns). The second function takesthe order argument indicating the order of the derivative to be taken.derivative – VISImage derivative(unsigned int order x, unsigned int order y) constVISImage derivative(unsigned int order x, unsigned int order y, float scale) constThese functions take the n-th order derivative in x and y directions. The derivatives can be of differentorder in the x vs. y direction. In the second function the scale parameter indicates the standarddeviation of a Gaussian kernel (see [5]) used to compute the derivative.derivatives – VISImage derivatives(unsigned degree) constVISImage derivatives(unsigned degree, float scale) constSimilar to the derivative functions above but return a multi-channel image with derivatives up to adegree.dxHalfForward/dxHalfBack/dyHalfForward/dyHalfBack – VISImage dxHalfForward() constVISImage dxHalfBack() constVISImage dyHalfForward() constVISImage dyHalfBack() constThese functions return “half” derivatives of an image, using forward or backward differences. Theseare useful in solving differential equations (see method anisoDiffuse()) These functions all assume atoroidal topology (periodic boundary conditions).zeroCrossings – VISImage zeroCrossings() constThis function places 0’s in pixel locations closest to zero crossings along grid lines. Other pixels in theimage are undisturbed.reduce – VISImage reduce(int scale) constReduces the size of an image by sub sampling, taking every scale’th pixel.reduceAverage – VISImage reduceAverage(int scale) constDoes the same as reduce but averages all pixels around the scale’th pixel and places the averagevalue in the returned image.resample – VISImage resample(float the scale) constVISImage resample(unsigned w, unsigned h) constVISImage resample(float scale x, float scale y, float x, float y) constThese functions scale the size of the image by interpolating between pixels. This change in scaling isspecified in three different ways. In the first function, the image size is changed by scale, which is afloating point value. In the second function, the image size is changed to match the size of w by h. Inthe third function the image size is scaled by scale x in width and scaled by scale y in height. Inthis function x and y indicate a position on the image upon which the center of the new image gridwill be placed.31noiseShot – VISImage noiseShot(float percent, T low, T high) constVISImage noiseShot(float percent, float sigma) constThis function adds shot noise to the image being operated on. The percent parameter represents thepercent of pixels effected by the noise. The low and high parameters are the limits of the uniformallydistributed noise and the sigma parameter is the sigma value of the Gaussian distributed noise.setBorder – VISImage setBorder(T value, unsigned thickness) constSets a strip of pixels of width thickness to value around the border of the image.becomeFlat – VISImage becomeFlat() constThis function creates a single channel image from a multichannel image by shrinking each channel andplacing the result into an image, that is, a ’matrix’ of the shrunken channel images. This ’matrix’image is the image returned.repeatDown – VISImage repeatDown(unsigned num) constThis function creates a bigger image by repeating copies of the image being operated on.repeatAcross – VISImage repeatAcross(unsigned num) constThis function creates a bigger image by repeating copies of the image being operated on.distanceTrans – VISImage distanceTrans() constVISImage distanceTrans(float max distance) constThese distance transforms return an image containing distance values. The distance values indicatethe Euclidean distance from the respective point in the input image to the closest non-zero point inthe input image. The second distance transform function has the input parameter max distance. Ifthe distance is greater than max distance then a zero is placed in the returned image.cityBlockDistTrans – VISImage cityBlockDistTrans() constCity block distance transform.extrema – VISImage extrema() constReturns an image of 1’s indicating local extrema and 0’s elsewhere.4.2.5 Image Processing AlgorithmsanisoDiffuse – VISImage anisoDiffuse(float k) constVISImage anisoDiffuse(VISImage image dx, VISImage image dy) constThese functions return one iteration of an anisotropic diffusion [11, 13]. The parameter k is the con-ductance parameter (expressed in units of RMS of the gradient image). Alternatively, one can computethe conductance in a separate function and pass the two x-offset and y-offset conductance functions,image dx and image dy respectively.watershed – VISImage watershed(float thresh, float depth, int pad)VISImage watershed(float thresh, float depth, int size, int pad)These functions apply the watershed algorithm to an image [2]. The thresh parameter is the valueof a threshold applied before performing the watershed. This parameter only serves to speed upthe algorithm. If set to zero then no pre-thresholding occurs. The depth parameter indicates themaximum depth of the watersheds created. The size parameter is the minimum size of a region.Regions containing more than size pixels are merged with other regions. The pad parameter is thesize of padding placed around the boarder of the image.floodFill – VISImage floodFill(T thresh, unsigned x, unsigned y)const VISImage& floodFill(T label from, T label to, unsigned x, unsigned y)The first function returns an image of ones indicating where the flood fill occurred. All other values32are unchanged. The x and y parameters are the starting “seed” position and the thresh parameterindicates the value of the maximum level filled. The second function does the flood fill and changesall values label from to label to (e.g., floodFill(1,5,4,4) means that all the 1’s connected to thepixel at 4,4 by a path of 1’s are replaced with 5’s).floodfileBoundBox – VISImage floodFillBoundBox(T thresh, unsigned x, unsigned y,unsigned &x lo, unsigned &y lo, unsigned &x hi, unsigned &y hi)This function is similar to the two floodFill functions, except that it returns a bounding box on theflooded region, through the parameters x lo, etc.cannyEdges – VISImage cannyEdges() constVISImage cannyEdges(T threshold) constThese functions return the zero crossings of the second derivative in the gradient direction (withnegative third derivative), i.e., non-maximal edge detection, closely related to the Canny edge detector[4]. The threshold parameter is a threshold applied to the gradient magnitude.5 Volume LibraryThe two main classes in the volume library (VL) are the VISVolume class and the VISVolumeFile class. TheVISVolumeFile class is used only for file I/O and the VISVolume class is used for all other functionality. Theheader files, volume.h and volumefile.h, which include these two classes, provide much inline informationabout the classes themselves. These header files in conjunction with this document should give users theinformation needed to create any type of application.Our discussion of the VL here consists of a short description of the functionality of the library followed bya reference section. The reference section is similar to the other reference sections in this document in thatit contains a listing of all the functions in the VL.Include files needed to use the library are:• volume.h• volumefile.h5.1 FunctionalityWe divide the functions of the VL in the same manner as the functions of the image library: general functionsand volume processing functions. The general functions include functions that allow the user to interactwith the data of a volume, file I/O functions, overloaded operators, and basic math functions. By volumeprocessing functions we mean functions that involve more than a basic operation. A volume processingfunction is built using the general functions of the VL. We break our discussion into these two parts.335.1.1 General FunctionsIn this section we cite sources of information on the general functions and provide a list describing thedifferent types of general functions of the VL 2. This information is available from more than one source.Detailed descriptions of the all general functions excluding functions involving file I/O exist found in thevolume.h file. The functions involving file I/O exist in volumefile.h. Section 5.2 of this document alsocontains descriptions of all the general functions.Accessing data – These functions involve ’getting’ and ’putting’ values from/to voxels. They use copy onwrite (section 2.1).File I/O – The VISVolumeFile provides methods for saving and loading volumes, including floating pointvolumes.Overloaded Operators – Many operators have been overloaded for the VISVolume class.Math Functions – Many math function that are not overloaded operators exist in the library.5.1.2 Volume Processing FunctionsGeneral volume processing functions exist that:• apply masks to a volume• smooth a volume• determine derivatives of a volume• resample a volume• perform distance transforms on a volume• corrupt a volume with noiseWe include detailed descriptions of these functions in the volume.h file. In addition to the functions describedin the above list, we have implemented the marching cubes algorithm. This algorithm is described in thevolumefile.h file. Descriptions of all of these functions and algorithms can be found in the next section.5.2 Reference for Volume LibraryThis reference section lists the functions of the VL. The functions are listed according to category.Table of contents:5.2.1 Interacting with the dataPeeking and PokingDirect Access to the Data2In addition section 4.1.2 discusses the general functionality of the image library. The general functionality of the imagelibrary is similar to that of the volume library, therefore it may be useful to look at section 4.1.2 if the list given in this sectionis not sufficient.34Constructing a VolumeCreating Images from VolumesFile I/O5.2.2 Overloaded Image Operators5.2.3 Math Functions5.2.4 General Volume Processing Functions5.2.5 Volume Processing Algorithms5.2.1 Interacting with the DataThis section includes several subsections entitled: Peeking and Poking, Direct Access to the Data, Con-structing a Volume, Creating Images from Volumes, and File I/O. All the functions listed here are membersof the VISVolume class unless otherwise noted.Peeking and Pokingpoke – T& poke(unsigned int x, unsigned int y, unsigned int z)Assign a value to the voxel of a volume indicated by x, y, and z. This function is always used on thelhs of an assignment. Example:VISVolume V(4,4,4);V.poke(2,3,3)=1.44;peek – const T& peek(unsigned int x, unsigned int y, unsigned int z) constReturns a value of the voxel of a volume indicated by x, y, and z. Example:VISVolume V(4,4,4);V.poke(2,3,3)=1.44;float m=V.peek(2,3);operator() – const T& operator()(unsigned int x, unsigned int y, unsigned int z) constSame as peek. Returns a value of the voxel of a volume indicated by x, y, and z. Example:VISVolume V(4,4,4);V.poke(2,3,3)=1.44;float m=V(2,3,3);interp – T interp(float x, float y, float z) constThis function is similar to peek and operator(), but it interpolates between voxel locations of thevolume to determine the value at the location represented by x, y, and z. This is a linear interpolation.getROI – VISImage getROI(unsigned int x pos, unsigned int y pos, unsigned int z roi,unsigned int w roi, unsigned int h roi, unsigned int d roi)This function returns a region of interest from a volume. The region has a position (x pos,y pos,z pos):its upper left corner. In addition the region has a size (w roi,h roi,d roi): its width, height, and depth.putROI – void putROI(const VISVolume& volume in, unsigned int x pos,unsigned int y pos, unsigned int z pos)This function allows a user to place a region of interest into a volume. The region itself is a volumeand has a position (x pos,y pos,z pos): its upper left corner.printData – void printData() constThis does a floating point printf on all of the data. Be careful if you have a large volume.35Direct Access to the DatarepRef – VISVolumeRep* repRef()This routine returns a non-const pointer to the data of a volumeConstructing a Volume Most of these functions are listed according to the input of the function.default – VISVolume()An empty volume is created.volume – VISVolume(const VISVolume& volume)This function creates a volume from another volume. The volume created is exactly the same asvolume.image – VISVolume(const VISImage& image)This function creates a volume from a multi-channel image. The volume created is exactly the sameas image.w,h,d – VISVolume(unsigned int w, unsigned int h, unsigned int d)This function creates a volume of a specified size, where w is width, h is height, and d is the depth,which correspond to the x, y, and z directions respectively.buf – VISVolume(unsigned int w, unsigned int h, unsigned int ch, T** buf)This function creates a volume w x h x d, using buf as the data. The buf buffer is a pointer to anarray of pointers to image-like buffers. This is used to create volumes from stacks of images.createToSize – VISVolume createToSize() constReturns a volume of the same size as the volume being operated on. Example:VISVolume V1(128,128,128);VISVolume V2;V2 = V1.createToSize; //V2 is now 128x128x128, values in the volume are unknownCreating Images from Volumesimage – VISImage image()VISImage image(unsigned int slice)The first function creates a multi-channel image that matches the volume. Each channel in one “slice”in depth (z direction). The second function creates a single-channel image corresponding to the numbergiven in slice.imageRef – VISImage imageRef()VISImage imageRef(unsigned int slice)The “ref” forms return images which contain actual pointers to the same buffers used by the volume.The resulting images are only valid as long as the associated volume, and the buffers are not cleanedup when the images are destroyed.File I/O The VISVolumeFile class contains all the file I/O functions including two read functions andtwo write functions. These functions read and write floating point and binary format files. In addition theclass includes a function that executes the marching cubes algorithm and saves the result in an ivview fileformat. We list this marching cubes function in our discussion of volume processing algorithms, rather thanhere, because this is not primarily a file I/O function.36read float and write float – VISVol read float(const char* fname)int write float(const VISVolume& volume, const char*) constThese functions read/write a ”raw” floating point file. This file has an ascii header. We use a ”raw”format because no good, widely-used volume file format exists. The header has the following format:volume.width() volume.height() volume.depth()Example (reading a floating point format file):VISVolume input; //input volumeVISImageFile file;if (!((input=file.read_float("volumein.vol")).isValid())){fprintf(stderr,"Error reading image file %s.\n","volumein.vol");exit(0);}input=VISVolume(input);//Cast input volume as a floating pointoutput = input + 5.0f;file.write_float(output,"volumeout.vol"); //write the volume as a "raw"//floating point fileread binary and write binary – VISVol read binary(const char*) constint write binary(const VISVolume& volume, const char*) constThese functions read/write a volume of type byte. The volume is assumed to be binary and the valuesare stored in the file as bits. The header is binary. The header has the following format(where thevalues are stored in the file as unsigned short’s):volume.width() volume.height() volume.depth()5.2.2 Overloaded Volume OperatorsIn some cases, functions may have the same name but different inputs and outputs. In these cases the labelsvolume and from indicate an object of type VISVolume; whereas, any other label ,such as value, refersto a value. Labels for values are type T unless otherwise noted.operator=(from) – VISVolume& operator=(const VISVolume& from)Assigns one volume to another. This function does not copy the data. (See section 2.1.) Example:VISVolume A(128,128,128);A=1.0f; //every pixel in A is 1VISVolume B;B=A; //B is the same as Aoperator=(value) – VISVolume& operator=(T value)Assigns each voxel in a volume to the value of value. Example:VISVolume A(128,128,128);A=3.141592;//Every element of A is equal to 3.141592operator+(volume) – VISVolume operator+(const VISVolume& volume) constAdds two volumes. The volumes must be the same size. Example:VISVolume A(128,128,128),B(128,128,128),C;A=1.0f;B=2.0f;C=A+B;//now C=3.0f37operator+(value) – VISVolume operator+(const VISVolume& from, T value)VISVolume operator+(T value,const VISVolume& from)Adds value to a volume.Example:VISVolume A(128,128,128),B;A=1.0f;B=A+2.0f;//now B=3.0fB=3.0f+A;//now B=4.0foperator-(volume) – VISVolume operator-(const VISVolume& volume) constSubtracts one volume from another volume. The volumes must be the same size. Example:VISVolume A(128,128,128),B(128,128,128),C;A=1.0f;B=2.0f;C=B-A;//now C=1.0foperator-(value) – VISVolume operator-(const VISVolume&from, T value)VISVolume operator-(T value,const VISVolume& from)Subtracts a value from a volume or subtracts a volume from a value. Example:VISVolume A(128,128,128),B;A=1.0f;B=2.0f-A;//now B=1.0fB=A-3.0f;//now B=-2.0foperator*(volume) – VISVolume operator*(const VISVolume& volume) constMultiplies each voxel in a volume to the respective voxel in volume. Example:VISVolume A(128,128,128),B(128,128,128),C;A=5.0f;B=2.0f;C=B*A;//now C=10.0foperator*(value) – VISVolume operator*(const VISVolume& from, T value)VISVolume operator*(T value, const VISVolume& from)Multiplies each voxel in a volume by value. Example:VISVolume A(128,128,128),B;A=2.0f;B=2.0f*A;//now B=4.0fB=A * -3.0f;//now B=-6.0foperator/(volume) – VISVolume operator/(const VISVolume& volume) constDivides each voxel in a volume by the respective voxel in volume. Example:VISVolume A(128,128,128),B(128,128,128),C;A=5.0f;B=2.0f;C=B/A;//now C=0.4operator/(value) – VISVolume operator/(const VISVolume& from, T value)VISVolume operator/(T value, const VISVolume& from)Divides each voxel in a volume by value or divides a value by each voxel in a volume. Example:VISVolume A(128,128,128),B;A=2.0f;B=10.0f/A;//now B=5.0fB=A/4.0f;//now B=0.538operator+=(volume) – VISVolume& operator+=(const VISVolume& volume)Adds each voxel in a volume by the respective voxel in volume and assigns the result to the volumebeing operated on.operator+=(value) – VISVolume& operator+=(const T& value)Adds value to each voxel and assigns the result to the volume being operated on.operator-=(volume) – VISVolume& operator-=(const VISVolume& volume)Subtracts each voxel in volume from the respective voxel in the volume being operated on and assignsthe result to the volume being operated on.operator-=(value) – VISVolume& operator-=(const T& value)Subtracts value from each voxel and assigns the result to the volume being operated on.operator*=(volume) – VISVolume& operator*=(const VISVolume& volume)Multiplies each voxel by the respective voxel in volume and assigns the result to the volume beingoperated on.operator*=(value) – VISVolume& operator*=(const T& value)Multiplies each voxel by value and assigns the result to the volume being operated on.operator/=(volume) – VISVolume& operator/=(const VISVolume& volume)Divides each voxel by the respective voxel in volume and assigns the result to the volume being operatedon.operator/=(value) – VISVolume& operator/=(const T& value)Divides each voxel by value and assigns the result to the volume being operated on.5.2.3 Math Functionspower(int power) – VISVolume power(const int& power) constThis function raises the volume to a power.sqrt() – VISVolume sqrt() constPerforms the square root of each voxel in a volume.div(volume, value) – VISVolume div(const VISVolume& volume, T valueA divide method that returns value for any voxel that is zero in volume.abs() – VISVolume abs() constReturns the absolute value of the volume.max() – T max() constReturns max value in volume.min() – T min() constReturns min value in volume.sum() – float sum() constReturns the sum of all voxels in the volume.average() – float average()Returns the average of all voxels in the volumescale(value) – VISVolume scale(float value) const;Multiplies each voxel by value.395.2.4 General Volume Processing Functionsgauss – VISVolume gauss(float sigma) constGaussian smoothing. The standard deviation of the Gaussian filter is sigma.convolve – VISVolume convolve(const VISVolume& kernel) constConvolution with any filter you choose.maskFloat – float maskFloat(const VISVolume& mask, unsigned int center x,unsigned int center y, unsigned int center z, unsigned int x pos,unsigned int y pos, unsigned int z pos) constArray* maskFloat(const Array< VISVolume >& masks,unsigned int x pos, unsigned int y pos const, unsigned int z pos) constDetermines the inner product of a mask and a region in a volume. In the first function, the center ofthe mask is indicated by center x, center y, center z and the center of the mask is positioned atx pos, y pos, z pos of the volume. The second function does the same as the first except it applies anarray of masks to the volume, and thus, returns an array of floats.dx – VISVolume dx() constVISVolume dx(unsigned int order) constDetermines the derivative in x direction (along width). The second function takes the order argumentindicating the order of the derivative to be taken.dy – VISVolume dy() constVISVolume dy(unsigned int order) constDetermines the derivative in y direction (along height). The second function takes the order argumentindicating the order of the derivative to be taken.dz – VISVolume dz() constVISVolume dz(unsigned int order) constDetermines the derivative in z direction (along depth). The second function takes the order argumentindicating the order of the derivative to be taken.derivative – VISVolume derivative(unsigned int order x, unsigned int order y,unsigned int order z) constThese functions take the n-th order derivative in x, y, and z directions. The derivatives can be ofdifferent order in the x, y, and z directions.zeroCrossings – VISVolume zeroCrossings() constThis function places 0’s in voxel locations closest to zero crossings along grid lines. Other voxels in theimage are undisturbed.reduce – VISVolume reduce(int scale) constReduces the size of a volume by sub sampling, taking every scale’th voxel.reduceAverage – Reduces the size of a volume by sub sampling, taking every scale’th voxel and averagingover a footprint of the same size (i.e. reduces aliasing artifacts as it subsamples).resample – VISVolume resample(float the scale) constVISVolume resample(unsigned w, unsigned h, unsigned d) constVISVolume resample(float scale x, float scale y, float scale z,float x, float y, float z) constThese functions scale the size of the volume by interpolating between voxels. The change in scaling isspecified in three different ways. In the first function, the volume size is changed by scale, which is afloating point value. In the second function, the volume size is changed to match the size of w x h x d.In the third function the volume size is scaled by scale x in width, scaled by scale y in height, and40by scale z in depth. In this function x, y, and z indicate a position on the volume upon which thecenter of the new volume grid will be placed.noise – VISVolume noise() constThis function returns a volume with random values drawn from a uniform distribution ranging zero toone.setBorder – void setBorder(T value)void setBorder(T value, unsigned w)The first function sets a strip of voxels (1 voxel wide) to value around the border of the volume. Thesecond function sets a strip of voxels of width w to value around the border of the volume.floodFill – VISImage floodFill(T thresh, unsigned x, unsigned y, unsigned z)This function returns a volume of ones indicating where the flood fill occurred. All other values areunaffected. The x, y, and z parameters are the starting “seed” position, and the thresh parameterindicates the value of the maximum level filled.cityBlockDistTrans – VISVolume cityBlockDistTrans() constCity block distance transform. Values returned are on the interval [0,1].extrema – VISVolume extrema() constReturns a volume of 1’s indicating local extrema and 0’s elsewhere.5.2.5 Volume Processing AlgorithmsWe have implemented only one volume processing algorithm. This function is a member of the VISVolumeFileclass because it creates a file as its output.write marchingCubes – write marchingCubes(const VISVolume& volume,float iso value, char* fname) constApplies the marching cubes algorithm to a volume [6]. The iso value parameter indicates the valuethat the algorithm uses to determine the surface location in the volume. The fname parameter isthe name of the ivview file created by the function. This file contains the surface resulting from themarching cubes algorithm.6 Util LibraryThe UL (util library) contains several classes that are simply utilities. These utilities support other parts ofVISPACK and can be used directly by users. In fact we encourage users to look through the functionalityof the UL before beginning to program with VISPACK in order to be aware of the functions available in theUL.6.1 ContentsThe UL contains two main classes and a collection of stand-alone functions. The two main classes are atemplated array class, Array, and a templated list class, List. The useful stand-alone functions inthe UL are found in th mathutil.h file. (In this section we refer to the functions in mathutil.h as mathfunctions.) In most cases the math functions in the UL have special functionality or special implementation(e.g., improved efficiency over traditional functions).416.2 Reference for Util LibraryThis reference section lists the functions of the UL. The functions are listed according to category.Table of contents:6.2.1Array ObjectConstructing ArraysDetermining SizeInterfacing with the Elements6.2.2List Object6.2.3 Math Functions6.2.1 Array ObjectHere we list the functions of the Array class. These functions are listed by category.Constructing Arrays Most of these functions are listed according to the input of the function.default – Array()An empty array is created.array – VISArray(const VISArray& s)This function creates an array from another array. The array created is exactly the same as array.size – VISArray(unsigned size0)This function creates an array containing a specified number of elements. The size0 parameter specifiesthe number of elements.x0, x1, x2, x3 – Array(T x0,T x1)Array(T x0,T x1,T x2)Array(T x0,T x1,T x2,T x3)These functions create a 2 element to 4 element sized array containing some of x0, x1, x2, and x3(depending on the function used).Determining Sizen – unsigned n() constThis function returns the number of elements in an array.size – unsigned size() constThis function does NOT return the number of elements in an array. We include this function in thereference only to point out that it should NOT be used to find the number of elements in an array.Example:Array A; //size()=0 n()=0A.appenditem(0.0f); //size()=8 n()=142Interacting with the Elementspeek – const T& peek(unsigned idx) constReturns the contents of element number idx in the array. (First element in the array is idx=0.)Example:Array m array(4);VISMatrix M(2,2);M=2.0f;m array.poke(2)=M;VISMatrix N;N=m array.peek(2);//Now N is the same as Moperator() – const T& operator()(unsigned idx) constSame as peek. Returns the contents of element number idx in the array. (First element in the arrayis idx=0.) For an example see peek.poke – T& poke(unsigned idx)Returns a reference to the data in element number idx of the array. This function is always used onthe lhs of an assignment.Example:Array m array(4);VISMatrix M(2,2);M=2.0f;m array.poke(2)=M;refAt – T& refAt(unsigned idx) constThis routine returns a const pointer to the data of a chosen element of an array.appendItem – void appendItem(const T& t)Adds an item onto the end of the array.prependItem – void prependItem(const T& t)Adds an item onto the beginning of the array.insertItemAt – void insertItemAt(const T& t,unsigned idx)Inserts an item at index idx.replaceItemWith – void replaceItemAtWith(unsigned idx, const T& t)Replaces an item in an array with a new item at index idx.removeItemAt – void removeItemAt(unsigned idx)Removes item from array at index idx.reverse – void reverse()void reverse(unsigned start,unsigned end)The first function reverses the order of the entire array. The second function reverses the elementsbetween indeces start and end.sort – void sort()This function performs a (descending) bubble sort, relying on the “¡” and “¿” operators for the templateclass object.operator= – Array& operator=(const Array& s)Assigns one array to another. This function copies the array.sizeTo – virtual void sizeTo(unsigned new size)Changes the size of an array to new size.43copy – Array* copy() constReturns a copy of an array.6.2.2 List ObjectThis list object is a simple implementation of a doubly-linked list. It consists of templated Link andList objects. This List is not a full-blown object, and must be handled by pointer (i.e. the copy isnot deep). The Link contains pointers to “before” and “after” objects on the list and a reference to data oftype T. The link has methods for access to data and inserting new elements either before or after, and hasa constructor which takes data of type T. The list contains a pointers to the head and tail and methods forinserting, removing, and getting access to specific links.The Link object contains the following methods:next – Link *next() constReturns a pointer to next link on list.prev – Link *prev() constReturns a pointer to previous link on list.append – virtual void append(Link* link)Inserts the link after the current one.prepend – virtual void prepend(Link* link)Inserts the link before the current one.data – virtual T& data()virtual void data(T data)Returns or sets data in link.The List object contains the following methods:head/tail Link *head() constLink *tail() constReturn the head/tail of the list, respectively.append/prepend Link *appendItem(T data)Link *prependItem(T data)Puts link containing data at the end/beginning of list, respectively.insertItem Link *insertItemBefore(T data,Link *link)Link *insertItemAfter(T data,Link *link)Puts a link containing data before/after the given link, respectively.6.2.3 Math Functionstabs – T tabs(const T& a)This function returns the absolute value of a.44power – float power(float a, int n)This function raises a to the power of n.min – const T& min(const T& a,const T& b)const T& min(const T& a,const T& b,const T& c)const T& min(const T& a,const T& b,const T& c,const T& d)const T& min(const T& a,const T& b,const T& c,const T& d,const T& e,const T& f)The first function returns a if ab else b is returned. The second function returns max(a,max(b,c)).The third function returns max(max(a,b),max(c,d)). The fourth function returns max(max(a,b,c,d),max(e,f)).rad to deg – VISReal rad to deg(VISReal radians)Converts radians to degrees.deg to rad – VISReal deg to rad(VISReal degrees)Converts degrees to radians.rand1 – float rand1()Returns a uniformly distributed deviate on the [0,1] interval.gasdev – float gasdev()Returns a normally distributed deviate with zero mean and unit variance, using ran1() as the sourceof uniform deviates.sqr – T sqr(T x)Returns the square of x (i.e. x2).cub – T cub(T x)Returns the cube of x (i.e. x3).gaussFast – float gaussFast(float a)Computes a normalized Gaussian on a using a look up table.gaussCumFast – float gaussCumFast(float a)Computes a normalized cumulative (integral of Gaussian) on a using a look up table.cosFast – float cosFast(float a)Fast cosine function using look up tables. Returns the cosine of a in radians.acosFast – float acosFast(float a)Fast arccosine function. Returns the arccosine of a in radians using a look-up-table.7 Level-Set Surface-Modeling LibraryThe Level-Set Surface-Modeling (LSSM) Library is an implementation of the level-set technique [8, 9] specif-ically for deforming surface models embedded in volumes. The implementation uses the sparse-field methoddescribed in [10]. The library implements all of the basic numerical algorithms and handles all of the datastructures required to perform LSSM. The strategy for using this library is to subclass the object VoxModel,45set some parameters, define a set of simple virtual functions that control the deformation process, initializethe model, tell the model to iteratively deform according to those equations. This document assumes thereader is somewhat familiar with level-set techniques, their numerical and LSSM as described in [12]. Webegin with a mathematical description of what we are doing and how it relates to the object in the library;this is a short review that does not substitute for a basic familiarity with the literature on the subject.We also present an example of using level-set models to do 3D shape metamorphosis as described in theliterature [10, 3].7.1 Surface DeformationThe LSSM library allows one to solve for surface deformations, as a function of time, for surface movementsof the following form:∂S∂t= αF (S,N(S)) + βG(S,N(S))N(S) + γN(S) + ηE (H(S),K(s)) , (2)where H(S) and K(S) are the mean and Gaussian curvatures, respectively, of the surface at a point S. Thisequation is solved by representing the surface as the Kth level set of an implicit function φ(x, t) : ℜ3×ℜ+ 7→ℜ. This gives∂φ∂t= αF (x,∇φ)) · ∇φ+ βG(x,∇φ)|∇φ|+ γ|∇φ|+ ηE(Dφ,D2φ), (3)where Dφ and D2φ are collections first and second derivatives of φ, respectively. In practice equation 3is solved on a discrete grid using an up-wind scheme gradient calculations, centralized differences for thecurvature, and forward finite differences in time. The updates for the nth iteration on the discrete grid uni,j,ktake the formun+1i,j,k = ∆uni,j,k∆t. (4)The LSSM library uses the sparse-field method [10], and therefore updates to not take place on the entiregrid ui,j,k at each iteration, but rather only on a finite subset of grid points that lie near the zero level set,which serves as the surface model. This set of grid points is called the active set.Thus, the LSSM library offers the following capabilities:1. Creates an initial model (with associated active set) from a volume.2. Calculates ∆uni,j,k and ∆t using virtual functions (defined by subclasses) that describe F and G, andparameters (values set by the subclass) α, β, γ, and η.3. Performs an update on the values of uni,j,k.4. Maintains the list of active grid points and updates the layers around those points in order to maintaina neighborhood from which to calculate subsequent updates.5. Provides access to the volume that defines uni,j,k and the linked list of active grid points.Given the volume defining uni,j,k, one can then rely on the functionality of the Volume Library (see Section 5)for subsequent processing, file I/O, or surface extraction.46Figure 2: The base class LevelSetModel keeps track of the active set and surrounding layers, while thederived class Voxmodel performs updates on the active set from a set of virtual functions.7.2 Structure and Philosophy of the LSSM LibraryThe library is organized (mostly for ease of development) into a base class, LevelSetModel, and a derivedclass, VoxModel, as shown in Figure 2. The base class does all of the book keeping associated with the activeset and surrounding layers, the link lists associated with those sets, and initializing the model. Thus it addsand removes voxels from the active set (and surrounding layers) in response to an update operation. Thebase class assumes that the subclasses know how to update individual voxels.The subclass, VoxModel, performs update on the grid points in the active set of the form given in Equation 4,using functions F and G and parameters α, β, γ, and η. It also calculates the maximum ∆t that ensuresstability.Thus, a user who wishes to perform a surface deformation using the LSSM library, would create subclassof VoxModel and define the appropriate virtual functions and set the parameters to achieve the desiredbehavior.Note, these two classes are not full fledged objects, but rather heavy-duty structure. The copyconstructors and assignment operators for these classes are not well behaved, and there-fore these classes should be passed by reference when sharing them between other rou-tines/modules.7.3 The LevelSetModel ObjectThe LevelSetModel contains a volume of values, a volume of status flags, five lists (one active list, two insidelists, and two outside lists), and three parameters that determine the origin of the coordinate system formwhich the model performs its calculations.There are two constructors, LevelSetModel() and LevelSetModel( const VISVolume &). Thefirst simply initializes the data structure, and the second also set the values of the model volume ( values)to the input. Once the values have been set, one can create an initial volume from those values by callingconstructLists(), which can also take a floating-point argument that controls the scaling of the inputrelative to a local distance transform near the zero set.The list that keeps track of the active set, called active list, keeps track of the location of those gridpoints and a single floating-point value, which stores the change in their values from one iteration to thenext.Another important methods for users of this object is update(float), which changes the grey-scale values ofthe grid for the active set according to the values stored in active list, and updates the status of elementson the active list as well as the values and status of nearby layers (2 inside and 2 outside). The floating pointargument is the value of ∆t from Equation 4, and the return value is the maximum change that occurredon the active set. Finally, the method iterate() calls the virtual method calculate change, a virtualfunction which sets the values of ∆uni,j,k and returns the maximum value of ∆t for stability, and then callsupdate. For this object the function calculate change performs some trivial (i.e., useless) operation.477.4 The VoxModel ObjectThe VoxModel object is a subclass of LevelSetModel, and it add three things to the base class.1. calculate change() is redefined to implement the surface deformation described in Equation 3.2. The virtual functions are declared for F (called force) and G (called grow). These functions aredefined to return null values for this object.3. Parameters are introduces which control the relative influence of the various terms. A routine load paramsis defined which simply reads these parameter from a file.4. A method rescale(float) is defined, which resamples the volume of grid-point values into a newvolume with different resolution and redefines the lists (and thereby the model) in this new volume.This method is for performing coarse-to-fine deformation procedures.7.5 Example: 3D Shape MetamorphasisThe Morph object allows one to construct a sequence of volumes or surface meshes using the 3D shapemetamorphasis technique described byBreen and Whitaker []. This technique relies distance transforms for both the source and target objectsand uses a LSSMs to manipulate the shape of the source so that it coincides with the target. The surfacedeformation that describes this behavior is∂S∂t= βG (T (S))N(S), (5)where G(x) is simply the distance transform (or some monotonic function thereof) of the target, and T is acoordinate transformation that aligns the source and target objects. The level-set formulation of this is∂φ(x, t)∂t= βG (T (S)) |∇φ|. (6)The morphing process consists of several steps:1. Read in distance transforms (in the form of volumes) for both source and target.2. Initialize the LSSM by fitting it to the zero set of the source distance transform.3. Update the LSSM according to Equation 6.4. Save intermediate volumes/surfaces at regular intervals.The remainder of this section lists the code and comments for three files, morph.h (which declares the Morphobject), morph.C (which defines the methods) and main.C (which performs all of the I/O and uses the Morphobject to construct a sequence of shapes.487.6 Morph.h//// morph.h////#ifndef iris_morph_h#define iris_morph_h#include "voxmodel/voxmodel.h"#include "matrix/matrix.h"#define INIT_STATE 0#define MORPH_STATE 1//// This is the morph object. It uses all of the machinery of the base// class to manipulate level sets. It needs to have an initial volume// and a final volume (which would typically be the distance transform,// it might need a 3D transformation, and it needs to redefine the// virtual function "grow", which takes 6 floats as input, the position// followed by the normal vectors (all will calculated and passed into// this method by the base class). It might also have a state, that// indicates whether or not it’s been initialized.//// Functions not defined here should be defined in "morph.C"//class Morph: public VoxModel{protected:VISVolume _dist_source;VISVolume _dist_target;VISMatrix _transform;//// This is the function that is used by the base class to manipulate thelevel// set. You can define it to by anything you want. For this object, itwill// return a value from the distance transform of the target.//virtual float grow(float x, float y, float z,float nx, float ny, float nz);// There are two states. In the first state, the model is trying to fit// to the input data. In this way the models starts by looking just like// the input dataint _state;public:Morph(const Morph& other)49{_dist_target = other._dist_target;_initial = other._initial;_state = MORPH_STATE;_transform = VISVISMatrix(3, 3);_transform.identity();// initialize();}Morph(VISVolume init, VISVolume d):VoxModel(){_dist_target = d;_initial = init;_state = MORPH_STATE;_transform = VISVISMatrix(3, 3);_transform.identity();// initialize();}void initialize();// for this object I assume that the transform is just a matrix.// but it could be anythingvoid transform(const VISVISMatrix& t){ _transform = t;}const VISVISMatrix& transform(){ return(_transform);}void distance(const VISVolume d){ _dist_target = d;}VISVolume distance(){ return(_dist_target);}};#endif7.7 Morph.C#include "morph.h"#include "util/geometry.h"#include "util/mathutil.h"//// this is the virtual function, that is the guts of it all.//float Morph::grow(float x, float y, float z,float nx, float ny, float nz)50{// this says you are in the morph state (things have been initialized)if (_state == MORPH_STATE){float xx, yy, zz;VISPoint p(4u);p.at(0) = x;p.at(1) = y;p.at(2) = z;p.at(3) = 1;VISPoint p_tmp;// this is where you could put some other transform.p_tmp = _transform*p;xx = p_tmp.x();yy = p_tmp.y();zz = p_tmp.z();// make sure you are not out of the bounds// of your distance volume.if (_dist_target.checkBounds(xx, yy, zz))// if not, get the distance (use trilinear interpolation).return(_dist_target.interp(xx, yy, zz));elsereturn(0.0f);}else{// if you are still initializing, then move toward the zero set of// your initial caseif (_initial.checkBounds(x, y, z))return(_initial.interp(x, y, z));elsereturn(0.0f);}}// this makes the model look like the input.#define INIT_ITERATIONS 5void Morph::initialize(){_values = _initial;int state_tmp = _state;_state = INIT_STATE;construct_lists(DIFFERENCE_FACTOR);// these couple of iterations are required to make sure that the zero// sets of the model match the zero sets of the//for (int i = 0; i < INIT_ITERATIONS; i++){// limit the dt to 1.0 so that the model settles in to a solution51update(::min(calculate_change(), 1.0f));}_state = state_tmp;}7.8 Main.C#include "vol/volume.h"#include "vol/volumefile.h"#include "image/imagefile.h"#include "morph.h"#include const int V_HEIGHT = (40);const int V_WIDTH = (40);const int V_DEPTH = (40);#define XY_RADIUS (12) // this matches the 2.5D data generated intorus.C#define T_RADIUS (4) // this matches the 2.5D data generated in torus.C#define S_RADIUS (12) // radius of a sphere#define B_WIDTH (20.0f)#define B_HEIGHT (60.0f)#define B_DEPTH (20.0f)#define B_CENTER_X (12.0f)#define B_CENTER_Y (32.0f)#define B_CENTER_Z (12.0f)float sphere(unsigned x, unsigned y, unsigned z);float torus(unsigned x, unsigned y, unsigned z);float cube(unsigned x, unsigned y, unsigned z);// This is a program that does the morph. If you give it two// arguments, it reads the initial model and the dist trans for the// final model from the two file names given, otherwise, it makes asphere// and deforms it into a torusmain(int argc, char** argv){VISVolume vol_source, vol_target;VISVolumeFile vol_file;int i;char fname[80];52vol_source = VISVolume(25,65,25);vol_source.evaluate(cube);if (argc > 2){// read in the sourceing modelvol_source = VISVolume(vol_file.read_float(argv[1]));// read in the dist trans of the final modelvol_target = VISVolume(vol_file.read_float(argv[2]));}else// make up some volumes{vol_source = VISVolume(V_WIDTH, V_HEIGHT, V_DEPTH);vol_source.evaluate(sphere);vol_target = VISVolume(V_WIDTH, V_HEIGHT, V_DEPTH);vol_target.evaluate(torus);}// create morph objectMorph morph(vol_source, vol_target);// loads in some parameters (for morphing these are all zero but one)// i.e.////////morph.load_parameters("morph_params");morph.initialize();vol_file.write_float(morph.values(), "morph0.flt");float dt;// do 150 iterations for your model to get from start to finish// probably don’t need this many iterationsfor (i = 0; i < 150; i++){dt = morph.calculate_change();// limit dt to 0.5 so that model never overshoots goaldt = min(dt, 0.5f);morph.update(dt);printf("iteration %d dt %f\n", i, dt);if (((i + 1)%10) == 0){// save every tenth volumesprintf(fname, "morph_out.%d.dat", i + 1);vol_file.write_float(morph.values(), fname);}}53// save a surface model (i.e. marching cubes).vol_file.march(0.0f, morph.values(), ‘‘morph_final.iv’’);printf("done\n");}54A Parameter File ObjectsA.1 IntroductionWhen implementing algorithms, it is often the case that a significant number of values are needed to ini-tialize various constants, specify input and output files, or otherwise tune the performance of the algorithm.Specifying increasingly more and more parameter values as command line arguments becomes unwieldy forboth the developer and the user of the software. The ParameterFile object is intended to simplify theprocess by defining a standard file format for specifying parameters and a simple protocol for reading thosefiles.This document discusses the ParameterFile object design and its user interface. Examples are provided.A.2 ImplementationThe ParameterFile object encapsulates the process of opening, reading, and parsing a file into the instan-tiation of a single high level object with a limited public interface. The implementation consists of threecomponents: a token scanner, a parser, and a syntax tree. The parser queries the token scanner for inputand constructs a syntax tree of objects whose root node is the ParameterFile object.The token scanner (yylex()) is created using the flex utility, a standard lexical analyzer generator. Theset of tokens T ,T = (identity, string, integer,decimal, rparen, lparen)recognized by the scanner are encoded in the file param.l as regular expressions.The parser is generated with the standard yacc utility and then compiled into a function yyparse(). Wheninvoked, yyparse() calls yylex() for input and attempts to identify syntactic structures of the grammarG3 = (V, T, P,ParameterFile) where V is the set of variables, V = (ParameterFile,Parameter,Value),and P is the set of productions:ParameterFile⇒ Parameter∗Parameter⇒ lparen identity Value∗ rparenValue⇒ ( string | integer | decimal )As productions in the grammar are identified, yyparse() constructs a “syntax tree” of objects representingthose productions. The tree is constructed by the parser from the bottom up. All the objects are subclassedfrom a pure virtual VPF::node object, thereby enabling run-time polymorphism of node object methods.A.3 Parameter file formatA parameter file as described by the grammar G, is an ASCII text file supplying a list of parameters of theform,3G is defined in the file param.y.55(parameter_name optional_value1 optional_value2 ... optional_valueN)The parameter naming conventions are similar to C language variable naming conventions. A parametername must begin with an alphanumeric character and may consist of any of the letters A-Z,a-z; numbers0-9; and the underscore character .The optional values may be either integer, decimal, or string values. Decimal numbers can be either floatingpoint or double precision4. A number value lacking a decimal point character is assumed to be an integervalue.A string value is specified by enclosing any characters (except the newline and " characters) in quotationmarks ".(parameter_name "string value number 1.1" "stri%$ng value number 2")The value types for a parameter can be mixed. For example,(mixed_values 45.23 1 0 2 3234.99999 "string")C++-style commenting is supported. C-style commenting is not supported. All characters after doubleforward slashes // up to a newline character are ignored.A.4 The ParameterFile objectIntroduction When implementing algorithms, it is often the case that a significant number of values areneeded to initialize various constants, specify input and output files, or otherwise tune the performance ofthe algorithm. Specifying increasingly more and more parameter values as command line arguments becomesunwieldy for both the developer and the user of the software. The ParameterFile object is intended tosimplify the process by defining a standard file format for specifying parameters and a simple protocol forreading those files.This document discusses the ParameterFile object design and its user interface. Examples are provided.Implementation The ParameterFile object encapsulates the process of opening, reading, and parsing afile into the instantiation of a single high level object with a limited public interface. The implementationconsists of three components: a token scanner, a parser, and a syntax tree. The parser queries the tokenscanner for input and constructs a syntax tree of objects whose root node is the ParameterFile object.The token scanner (yylex()) is created using the flex utility, a standard lexical analyzer generator. Theset of tokens T ,T = (identity, string, integer,decimal, rparen, lparen)recognized by the scanner are encoded in the file param.l as regular expressions.4Most variants on these types are also supported (long double, etc). Check the header files for a complete list of typessupported.56The parser is generated with the standard yacc utility and then compiled into a function yyparse(). Wheninvoked, yyparse() calls yylex() for input and attempts to identify syntactic structures of the grammarG5 = (V, T, P,ParameterFile) where V is the set of variables, V = (ParameterFile,Parameter,Value),and P is the set of productions:ParameterFile⇒ Parameter∗Parameter⇒ lparen identity Value∗ rparenValue⇒ ( string | integer | decimal )As productions in the grammar are identified, yyparse() constructs a “syntax tree” of objects representingthose productions. The tree is constructed by the parser from the bottom up. All the objects are subclassedfrom a pure virtual VPF::node object, thereby enabling run-time polymorphism of node object methods.Parameter file format A parameter file as described by the grammar G, is an ASCII text file supplyinga list of parameters of the form,(parameter_name optional_value1 optional_value2 ... optional_valueN)The parameter naming conventions are similar to C language variable naming conventions. A parametername must begin with an alphanumeric character and may consist of any of the letters A-Z,a-z; numbers0-9; and the underscore character .The optional values may be either integer, decimal, or string values. Decimal numbers can be either floatingpoint or double precision6. A number value lacking a decimal point character is assumed to be an integervalue.A string value is specified by enclosing any characters (except the newline and " characters) in quotationmarks ".(parameter_name "string value number 1.1" "stri%$ng value number 2")The value types for a parameter can be mixed. For example,(mixed_values 45.23 1 0 2 3234.99999 "string")C++-style commenting is supported. C-style commenting is not supported. All characters after doubleforward slashes // up to a newline character are ignored.The ParameterFile object The ParameterFile object is a “nice” object, meaning that the copy con-structor and assignment operators are defined to perform deep copies. Instances of ParameterFile objectscan be safely passed by reference or by value.A ParameterFile object can be instantiated with a character string, a copy constructor, or by assignment,5G is defined in the file param.y.6Most variants on these types are also supported (long double, etc). Check the header files for a complete list of typessupported.57VPF::ParameterFile p1("filename"); // filenameVPF::ParameterFile p2(p1); // copy constructorVPF::ParameterFile p3; p3 = p2; // assignmentNote the use of the scope resolution operator. The ParameterFile objects are enclosed in a namespaceVISParameterFile aliased to VPF.When instantiated with a character string argument, the constructor attempts to open and parse the filewhose name is represented in the argument, calling the yyparse() method. On success, a pointer to theresulting syntax tree is held in the ParameterFile object. The ParameterFile object contains the followingpublic methods:bool valid() Returns true if object contains a valid syntax tree,false otherwise.int size() Returns the number of parameters contained in the syntax tree.bool empty() Returns true if the object’s valid syntax tree is empty,false otherwise.void clear() Clears the syntax tree and frees all associated memory.void print() Prints a formatted list of everything in the syntax tree tostandard out. For debugging.The syntax tree of the ParameterFile object can be thought of as an ordered list of Parameter objects. The[] operator has been overloaded to allow (read) access to the parameters. Supplying an integer argumenti returns the ith parameter in the list, as ordered in the input file. No bounds checking is performed bythis operation. Supplying a character string argument returns the the parameter whose name matches theargument. If no match is found, a “null parameter” is returned whose validity is false (see next section).The Parameter object A Parameter object is also a “nice” object in the sense that it can be passed byvalue or reference and its assignment operator is overloaded. A Parameter object has the following methods:bool valid() Returns true if the parameter is valid, false otherwise.int size() Returns the number of values contained in the parameter.bool empty() Returns true if the object’s value list is empty,false otherwise.const char *getName() Returns the name of the parameter.void print() Prints the parameter name and its values tostandard out. For debugging.A Parameter object contains a character string name and an ordered list of Value objects. The [] operatorhas been overloaded to allow (read) access to this list. No bounds checking is performed by the [] operator.Note that since a Value object may be of several types, the assigment operator cannot be overloaded forthe Value object. Hence the value contained in the object can only be accessed by downcasting to theappropriate class and using that class’s getValue method (see example section). The preferred method foraccessing a value is through the overloaded set method described in the next section.The set method The VPF::set method has been overloaded for all the supported variable types. It takesas its arguments a variable of the appropriate type and a Value object and returns true if the operation58// This is a sample file to demonstrate the format readable by the// ParamterFile object// Parameters are enclosed in parentheses and take the form://// (parameter_name optional_value1 optional_value2 .... optional_valueN)//// C++ style comments can be used. Note C-style comments are NOT supported.//(bounding_box1 0 100 50 -50) (k .15)(forsythia "the tree with the" "lights in it" "edgy" "2394 sljdkf09128 ...sdfk238""32ls lskdl 02302...&*&%$#!" "right now strings cannot be more""than VPF::MAXSTRLEN characters long" )(flag_value) // No parameter values is OK, too.(mixed_types 0 "string_value" 34.2222 -0.987) // mixes types are OK, too// SOME COMMENTS(d 1.13425835) // MORE COMMENTSFigure 3: sample.txtwas successful and false otherwise. The declaration of the set method has the following form:VPF::bool set( &, VPF::Value &)The method works by first checking the type of the supplied Value object to make sure it is consistant withthe requested variable type. Then the value contained in the Value object is cast to the appropriate typeand the assignement is made.Error handling At the time of this writing, C++ exception handling is not used in the ParameterFilecode. Some run-time warnings are printed to standard error and true exceptions such as “File not found”result in a call to exit(). Syntax and scanner errors are reported to standard error, but the parser cannotusually recover from such errors and will cause a ParameterFile instantiation to fail with a call to exit.Examples Figure 3 shows a sample input file that is readable by the ParameterFile parser. Figures 4shows a small sample program that reads in some of these values.B Code ExamplesThe code examples shown here can be found in the examples directory. There are three examples: matrixuses the Matrix Library, shotnoise uses only the Image Library, and fractal uses both the Volume andImage Libraries. Examine the README files in the directories of these examples to learn the function of59#include #include "param.h"main(int argc, char *argv[]) {int x1, x2, y1, y2;float k;double d;char string[VPF::MAXSTRLEN];VPF::ParameterFile F("sample.txt");if ( ! F.valid() ) { cerr << "Problems reading input file." << endl; exit(1); }if ( F.empty() ) { cerr << "Input file is empty" << endl; exit(2); }if (F["bounding_box1"].valid() && F["bounding_box1"].size() == 4){ if ( ! VPF::set(x1, F["bounding_box1"][0]) )cerr << "Argument 0 of \"bounding box\" is not an integer value" < 0){ if ( !VPF::set(k, F["k"][0]) )cerr << "Argument 1 of \"k\" is not a decimal value" << endl;}else { cerr << "\"k\" is not a valid parameter name" << endl; }// Set a double valueif (F["d"].valid() && F["d"].size() > 0){ if ( !VPF::set(d, F["d"][0]) )cerr << "Argument 1 of \"d\" is not a decimal value" << endl;}else { cerr << "\"d\" is not a valid parameter name" << endl; }// Set a string value.// Note that the string set function assumes a buffer of at least VPF::MAXSTRLEN// characters.if (F["forsythia"].valid() && F["forsythia"].size() > 4){if ( !VPF::set(string, F["forsythia"][3]) )cerr << "Argument 3 of \"forsythia\" is not a string value" <#include #include "matrix.h"/////////////////////////////////////////////////////////////////////////// Summary: This program contains several ’main’ functions. Each //// one tests a different aspect of the matrix library. //// By commenting/uncommmenting portions of the code, //// different aspects of the matrix library may be examined //// Output: usually stdout// Author: Samuel G. Burgiss Jr. /////////////////////////////////////////////////////////////////////////////test of operator=main(char** argv, int argc){cout << "m1" << endl;VISMatrix M1(3,3);M1=1.0f;cout << "m2" << endl;VISMatrix M2;cout << "m2made" << endl;M2=M1;cout << "m3" << endl;VISMatrix M3;M3=3.0f;cout << "m3made" << endl;M1=M3;cout << M2 << endl;}//test of det()/*main(char** argv, int argc){VISMatrix A("input4");cout << A << endl;cout << A.det() << endl;61}*///test of operator()/*main(char** argv,int argc){VISMatrix A("input4");cout << A(1,1) << endl;cout << A(2,2) << endl;}*///test of operator==()/*main(char** argv,int argc){VISMatrix A("input4");VISMatrix B,C,D;B = A;B.poke(1,1) = 100.0f;C = A==A;D = A==B;cout << A << endl;cout << B << endl;cout << C << endl;cout << D << endl;}*//*//test pokeROI and peekROImain(char** argv,int argc){VISMatrix A("input6");VISMatrix B("input5");cout << "A" << endl << A << endl;cout << "B" << endl << B << endl;cout << "peekROI(0,2,4,5) - startrow,endrow,startcol,endcol" << endl;cout << A.peekROI(0,2,4,5) << endl;A.pokeROI(1,1,B);cout << "A.pokeROI(1,1,B);" << A << endl;VISVector a;a=A.vec(2);cout << "a" << endl << a << endl;cout << "a.peekROI(2,3)" << endl << a.peekROI(2,3) << endl;VISVector b;b = a.peekROI(2,3);cout << a.pokeROI(8,b) << endl;}*///test of +=, -=, *=/*main(char** argv,int argc){VISMatrix A("input4");VISMatrix B("input5");cout << "A" << endl << A << endl;cout << "B" << endl << B << endl;A+=B;cout << "A+=B" << endl << A << endl;cout << "B" << endl << B << endl;62A+=A;cout << "A+=A" << endl << A << endl;A-=A;cout << "A-=A" << endl << A << endl;A+=B;A+=B;cout << "A+-B;A+=B" << endl << A << endl;A*=A;cout << "A*=A" << endl << A << endl;cout << "A/5.0f" << endl << (A/5.0f) << endl;A/=5.0f;cout << "A/=5.0f" << endl << A << endl;A*=5.0f;cout << "A*=5.0f" << endl << A << endl;}*///test for vector class/*main(char** argv, int argc){VISVector v1;VISMatrix m(6,1);for (int i=0;i<6;i++)m.poke(i,0)=(float)i;cout << "Assignment and creation of a vector." << endl;cout << "m:" << endl << m << endl;v1=m+m;cout << "v1=m+m:" << endl << v1 << endl;cout << "Dot product." << endl;cout << "v1 dot v1: " << v1.dot(v1) << endl;VISVector v2(3);cout << "v2" << endl << v2 << endl;VISVector v3(3);cout << "crap" << endl;v2=1.0f;cout << "v2" << endl << v2 << endl;v2.poke(1)=2.0f;cout << "v2:" << endl << v2 << endl;v3=3.0f;cout << "v3:" << endl << v3 << endl;v2.cross(v3);cout << "v2 X v3:" << endl << v2.cross(v3) << endl;63cout << v2 << v3 << endl;VISVector v4;cout << "v2+v3" << endl << v2+v3 << endl;v4 = v2+v3;cout << "v2" << endl << v2 << endl;cout << "v3" << endl << v3 << endl;cout << "v2+v3" << v2+v3 << endl;cout << "ans" << endl << v4 << endl;}*///test the concat algorithms/*main(char** argv, int argc){VISMatrix a("input");cout << "a" << endl << a << endl;VISVector b(4);b=2.0f;b.poke(0) = 3.0f;b.poke(1) = 4.0f;cout << "b" << endl << b << endl;cout << "a.concatcol(b)" << endl << a.concatcol(b) << endl;cout << "b.concatrow(b)" << endl << b.concatrow(b) << endl;cout << "a.concatrow(a)" << endl << a.concatrow(a) << endl;}*///testing of the svd decomposition and//the covariance function//fitting a plane/*main(char** argv, int argc){VISMatrix C,U,W,V;VISVector n;float minval,d;int r;VISMatrix b("pts3");C=(b.t()).cov();cout << "b.cov()" << endl << C << endl;C.svd(U,W,V);cout << "U" << endl << U << endl;cout << "W" << endl << W << endl;cout << "V" << endl << V << endl;minval = (W.vecfromdiag()).min(r);cout << "minval " << minval << endl;cout << "r " << r << endl;n = -(V.vec(r));cout << "n" << endl << n << endl;d = -1.0f*((n.t()*(b.t()).mean()).peek(0));cout << "d " << d << endl;}*/64//test the min,max stuff./*main(char** argv, int argc){VISMatrix b("input3");cout << "b" << endl << b << endl;int r,c,n;cout << "b.max() " << b.max() << endl;n=b.max(r,c);cout << "n " << n << endl;cout << "r " << r << endl;cout << "c " << c << endl;}*///test transpose/*main(char** argv, int argc){VISVector v(1,2,3);cout << "v" << endl << v << endl;cout << "v.t()" << endl << v.t() << endl;}*/// testing the fstream stuff// writing to and reading from a file/*main(char** argv, int argc){VISMatrix temp_matrix;ifstream ultrafile("input");if(!ultrafile){cout << "Unable to open the file: input" << endl;exit(0);}ultrafile >> temp_matrix;cout << "temp_matrix" << endl << temp_matrix << endl;}*///test inverse functions/*main(char** argv, int argc){float **nra,**nrb;int i,j;int m=3;int n=3;VISMatrix a(m,n),b;a.poke(0,0)=2.0f;a.poke(0,1)=4.0f;a.poke(0,2)=2.0f;a.poke(1,0)=13.0f;a.poke(1,1)=23.0f;a.poke(1,2)=24.0f;a.poke(2,0)=4.0f;a.poke(2,1)=23.0f;a.poke(2,2)=7.0f;65cout << "a: " << endl << a << endl;b = a.inversesvd();cout << "b: " << endl << b << endl;cout << "inversegj" << endl << a.inversegj() << endl;}*/\end{tt}B.2 shotnoise example\begin{tt}#include "image/image.h"#include "image/imagefile.h"#include "util/mathutil.h"#include #include #include #include /////////////////////////////////////////////////////////////////////////// Summary: Adds shot noise to an image. //// Usage: shotnoise percent min max input_filename output_filename//// Output: fits file (filename chosen by user) //// Author: Samuel G. Burgiss Jr. ///////////////////////////////////////////////////////////////////////////int process(int num);main(int argc, char** argv){char infilename[80];char outfilename[20];float min,max;FILE *fp;VISIm image;VISImage output_image;VISImage input_image,deriv_x,deriv_y,added;VISImageFile input_file,output_file;float percent;if (argc<4){cout << "summary: Adds shot noise of particular characteristics" << endl;cout << " to an image." << endl;cout << "usage: shotnoise percent min max input_filename output_filename" << endl;cout << " percent = percent of pixels affected" << endl;cout << " min = min of the noise amplitude" << endl;cout << " max = max of the noise amplitude" << endl;exit(0);66}// Load in the argumentspercent = (float)atof(argv[1]);min = (float)atof(argv[2]);max = (float)atof(argv[3]);strcpy (infilename,argv[4]);strcpy (outfilename,argv[5]);// Open the file and check to see that it is valid.fp = fopen(infilename, "r");if (fp == NULL){cout << "Unable to open " << infilename << endl;exit(0);}if (!((image=input_file.read(infilename)).isValid())){cout << "Error reading image file " << infilename << endl;exit(0);}input_image=VISImage(image);// Print out the characteristics of the noise.cout << "percent = " << percent << endl;cout << "min = " << min << endl;cout << "max = " << max << endl;// Add the noise to the imageoutput_image=input_image.noiseShot(percent, min, max);// Write the output file in fits formatoutput_file.write_fits(output_image,outfilename);}\end{tt}B.3 fractal example\begin{tt}#include "image/image.h"#include "image/imagefile.h"#include "volume.h"#include "volumefile.h"#include #include #include #include /////////////////////////////////////////////////////////////////////////// Summary: Creates a 3D surface by applying the marching cubes //// algorithm to a volume of data that is created using the //// mandelbroth set. //// Output: fits file (filename chosen by user) //// Author: Samuel G. Burgiss Jr. ///////////////////////////////////////////////////////////////////////////67void write_volume_flat(VISVolume V, char* filename);main(int argc, char** argv){float tempza;float begin_c,begin_ci;float end_c,end_ci;int number_slices;float begin_w;float end_w;int window_w;float begin_h;float end_h;int window_h;float increment_w,increment_h,increment_c,increment_ci;int volume_index_w,volume_index_h,slice_index;float cnt_w,cnt_h,cnt_c,cnt_ci;float za,zb;int iter,max_iter;float march_thresh;begin_c = -2.0f;begin_ci = 0.0f;end_c = 2.0f;end_ci = 0.0;number_slices = 100;begin_w = -2.0f;end_w = 2.0f;window_w = 100;begin_h = -2.0f;end_h = 2.0f;window_h = 100;max_iter = 30;march_thresh = 5.5;VISVolume V(window_w+3,window_h+3,number_slices+3);V=0.0f;increment_w = (end_w-begin_w)/(float)window_w;increment_h = (end_h-begin_h)/(float)window_h;increment_c = (end_c-begin_c)/(float)number_slices;increment_ci= (end_ci-begin_ci)/(float)number_slices;volume_index_w=1;slice_index=1;cnt_ci=begin_ci;cout << "Making fractal volume";for (cnt_c=begin_c;cnt_c V, char* filename){int slice,row,col;int v_width,v_height,v_depth;VISImage image_volume(V.width(),V.height(),V.depth());VISImage temp_image(V.width(),V.height());VISImageFile im_file;v_width = V.width();v_height = V.height();v_depth = V.depth();for(slice=0;slice

Vispack

Abstract

Similar works

Full text

Available Versions

The University of Utah: J. Willard Marriott Digital Library