Report

Introduction to MATLAB and MATLAB Programming #1 Outline • Who are we and why are we here? (to paraphrase James Stockdale, 1992) • Who am I and how did I get here? • Who are you and why are you here? • MATLAB • Started programming (BASIC) in 1981 • Started imaging algorithm and software development (Pascal) in 1983 • Cambridge Instruments Quantimet 900 – 1 MB memory – Booted from 8” floppy disk – $100,000 • BCC ASET (1983), AA (1986) – BASIC, Pascal, FORTRAN, x86 Assembler • TSU BSEE (1986) • NCSU MSEE (1988) – MATLAB – C/C++ – APL, ProLog • NCSU PhD (1996) 3D ICUS • 3D Intracoronary Ultrasound • Lumen and medial-adventitial boundary estimation 3D Imaging Geometry xw y ~ X w zw 1 zw xw yw v u u ~ x v 1 C u v u ~ x v 1 C TIPS Surgical Planning MRA Vessel Extraction Objectives • Develop a basic working knowledge of MATLAB (notation, syntax, built-in functions) • Learn to use the command line (MATLAB as a calculator) • Learn to develop scripts and functions • Learn basic programming skills • Learn to use the debugger to troubleshoot scripts and functions Objectives • Generate 2D and 3D data plots • Read, write, and manipulate numeric data (e.g. text files, binary files, Excel) • Read, write, manipulate, and display image data • Learn the basics of linear algebra • Learn about the FFT and complex variables What is MATLAB? • MATLAB - Matrix Laboratory • What is a matrix? – Fundamental data type in MATLAB – An image • What is an image? – A rectangular array of numbers – A point in a high-dimensional space – A representation of a system of equations What is MATLAB? • What is a matrix? n m n m – A linear mapping from R to R (or C to C )** • What’s a linear mapping? • Lines map to lines. Origin maps to the origin. – A collection of vectors – What is a vector? » Magnitude and direction (velocity, force, etc.) » An n-tuple of numbers What is MATLAB? • History – Netlib (netlib.org) • Large software library written in FORTRAN (FORmula TRANslator) – Cleve Moler (EISPACK, LINPACK, cofounder of MathWorks) • Computing environment for technical computation, visualization, design, simulation, and implementation in a wide range of application areas. Common MATLAB Applications • Numerical Methods – Linear algebraic equations – Roots and optimization – Curve fitting – Integration and differentiation – Differential equations • MATLAB toolboxes – Statistics, optimization, image processing, computer vision, bioinformatics, … Why MATLAB? • Pros: – Easy to start using MATLAB with little or no prior programming experience – Relatively easy to learn – Can be used as a calculator or as a programming language – Can be used to test simple algorithms (command line or scripts) prior to coding in a higher-level language Why MATLAB? • Pros: – MATLAB has a large library of optimized, robust mathematical functions – MATLAB has a large collection of very powerful toolboxes (e.g. image processing, statistics) – Easy to produce 2D and 3D graphs Why Not MATLAB? • Cons: – Cost associated with ease of use • Interpreted language • Slow (10-100x) compared to some other programming languages (e.g. C, C++, FORTRAN, Pascal) • May be too slow to solve large problems – Fairly expensive (free alternatives available) – Designed for scientific computing • Not always the best solution for other applications 2-by-2 System of Equations • Can solve graphically • Can add/subtract equations x y 4 x y 0 2x 4 2y 4 x2 y 2 2-by-2 System of Equations • Can solve with the cross or outer product – cross () 1 1 4 2 2 1 1 4 2 2 4 0 2 1 n-by-n System of Equations • What about higher dimensions? • Cross product only works in 3D. • Harder (impossible?) to visualize. – Circuit analysis example • Graph theory (non-planar) and linear algebra • 6-by-6 system of equations Matrix Operations • The entry (or element) in row i and row j of a matrix, A, is denoted by aij • The entry in row i of a vector (i.e. matrix with only one column) v is denoted by vi. a11 A a 21 a 31 a12 a13 a 22 a 23 a 32 a 33 a14 a 24 a 34 v1 v2 v v3 v4 Matrix Operations • Two matrices are equal if: – They are the same size – All corresponding entries are equal • Two matrices can be added together if they have the same size (i.e. the same number of rows and columns) C = A+B c ij a ij bij Matrix Operations • Multiplication by a scalar a11 cA c a 21 a12 ca11 a 22 ca 21 ca12 ca 22 • Matrix-vector multiplication n b Ax bi a ij x j j 1 – If A is m-by-n and x is n-by-1, b is m-by-1. Matrix Operations • Matrix-vector multiplication Ax = b a11 a 21 a 31 a12 a13 a 22 a 23 a 32 a 33 x1 a14 a11 x1 a12 x 2 a13 x 3 a14 x 4 b1 x2 a 24 a 21 x1 a 22 x 2 a 23 x 3 a 24 x 4 b 2 x 3 a 34 a 31 x1 a 32 x 2 a 33 x 3 a 34 x 4 b3 x4 Matrix Operations • Matrix-matrix multiplication n C AB c ij a ik b kj k 1 – If A is m-by-r and B is r-by-n, C is m-by-n. Matrix Operations • Matrix-matrix multiplication – To find cij (i.e. the element of C in row i and column j) multiply corresponding entries in row i of A and column j of B, and add them together. C AB a11 a 21 a12 b11 a 22 b 21 b12 a11b11 a12 b 21 b 22 a 21b12 a 22 b 21 a11b12 a12 b 22 c11 a 21b12 a 22 b 22 c 21 c12 c 22 Matrix Operations AB BA A + B + C = A + B + C A BC = AB C A B + C = AB + AC B + C A = BA + CA AB BA Matrix Operations • Zero matrix – A matrix where each element is 0 – zeros() 0 0 04 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Matrix Operations • Identity matrix – Square (i.e. m=n) – Diagonal entries equal to 1; all others are 0. – Multiplication by the identity matrix doesn’t change the matrix – eye() 1 0 0 0 0 I4 0 0 1 0 0 1 0 0 0 0 1 Row Reduction • Elementary row operations – Multiply a row by a scalar. – Add a scalar multiple of one row to another. – Exchange two rows. • Solve Ax=b – Row reduction with augmented matrix a11 a 21 a 31 a12 a13 a 22 a 23 a 32 a 33 b1 1 b2 0 0 b3 0 0 1 0 0 1 x1 x2 x 3 2-by-2 cont’d • Back to 2-by-2 example – Rewrite as a matrix equation x y 4 x y 0 1 1 1 x 4 1 y 0 – Form an augmented matrix 1 1 1 1 4 0 2-by-2 cont’d • Solve by row reduction – Make left side look like an identity matrix using only elementary row operations 1 1 1 0 1 1 1 2 4 0 4 4 1 1 0 1 4 2 1 0 2 2 0 1 Rank • r=rank(A) • Number of leading 1’s in the row-reduced matrix • From previous example, rank(A)=2 1 1 1 1 4 1 0 0 0 1 2 2 Matrix Inverse • Matrix inverse – A must be square and full rank – If a matrix can be found such that AB=In, then B is called the inverse of A. 2-by-2 cont’d • Computing the inverse A 1 1 1 0 I I 1 1 1 0 1 0 0 1 1 1 2 1 1 1 0 1 -1 A 1 0.5 0 0.5 1 0.5 0 1 0 0.5 0.5 0.5 1 x A b 0.5 0.5 0.5 4 2 0.5 0 2 Matrix Inverse -1 • If A is n-by-n and rank(A)=n, then A exists. – Implies there is a solution for every b • b is in the range of A – inv() • If A and B are both invertible, 1 A 1 AB A 1 1 B A 1 2-by-2 cont’d • Inconsistent (no solution) – n=2, Rank=1 – Inverse doesn’t exist x y 4 x y 5 1 1 1 0 1 x 4 1 y 5 1 0 4 1 2-by-2 cont’d • Consistent (infinite solutions) – n=2, Rank=1 – Inverse doesn’t exist x y 4 2x 2 y 8 1 2 1 x 4 2 y 8 1 0 1 0 4 0 Matrix Transpose • Transpose (‘ in MATLAB) – Swap rows and columns B A T A T bij a ji 1 3 5 A A B A + B T cA T 2 1 4 2 6 T 3 4 5 6 T AB T T cA T B A T T T Scalar Product • Also called inner or dot product – Vectors must have the same length (size) – Multiply corresponding elements and add n uv = uv i i 1 i u 1 v1 u 2 v 2 u nvn Scalar Product • Inner (dot) product u1 u = u2 u 3 v1 v v2 v 3 u v u 1 v1 u 2 v 2 u 3 v 3 T – Also written as: uv u v u T v cos Scalar Product • Properties u v v u u v w u v u w k u v ku v u kv vv 0 vv 0 v = 0 Norms • Vector norms – 1-norm n x 1 xi i 1 – 2-norm (Euclidean norm or length) n x 2 xi 2 T x x i 1 – Infinity-norm x m ax x1 , x 2 , , xn Norms • Matrix Norm – 1-norm (column sum) m A 1 m ax a ij 1 j n i 1 – 2-norm – later n – Infinity norm (row sum) – Frobenius norm A m ax a ij 1 i m m A F j 1 n i 1 j 1 2 a ij Norms • Vector norm properties x 0 x 0 x0 x x x+y x y • Matrix norm properties A 0 A A A+B A B AB A B Vector Product • Also called the cross product – 3D only u 2 v 3 - u 3 v 2 u v u 3 v1 - u1v 3 u 1 v 2 - u 2 v 1 Vector Product • Properties u v u v uv w u v u w u v w u w v w k u v ku v u kv u0 0u 0 uu 0 Unit Vectors • Unit vector v uˆ v uˆ 1 • Projections T w u v v 2 v Noise • Always present • Quantization noise • Finite word length (data type differences (intn, uintn, float, double) – Catastrophic cancellation • Electronic noise, photon noise, etc. • When is it a problem? Depends on the problem. • Be aware of effects of noise! Ill-Conditioning x 2y 4 2 x 3.999 y 7.999 1 2 x 4 3.999 y 7.999 2 Ill-Conditioning • Results are sensitive to noise • Loss of precision • cond() – condition number x 2 y 4.001 2 x 3.999 y 7.998 1.001 x 2.001 y 4 2.001 x 3.998 y 7.999 Ill-Conditioning • Trefethen least-squares ill-conditioned example (14th-order polynomial fit) • Solution by different methods – \ - matrix divide – inv() – matrix inverse – LU() – LU decomposition – QR() – QR decomposition – SVD() – Singular value decomposition MATLAB Workspace • • • • • Using MATLAB as a calculator Command window Command history Current directory Workspace MATLAB Workspace • • • • Help browser diary() Saving/loading the workspace Quitting MATLAB Variables, Assignments, and Keywords • Datatypes - scalars, vectors ([1, 2, 3, 4], [1:4]), matrices ([m,n], [1, 2; 3, 4]), uintn, intn, double, char, structs, cells, classes • Vector – magnitude and direction, n-tuple n m n m) • Matrix – mapping from R to R (or C to C , collection of vectors, representation of an image, representation of a system of equations Variables, Assignments, and Keywords • • • • Entering variables, vectors and matrices Incrementation and overwriting of variables Recalling expressions and making corrections Addressing vector and matrix elements – v(m:n), v(m:step:n), v([1 3 5 2 4 6]) – A(:,j), A(i,:), A(i,j:end) Matrix Properties • • • • • • • • • size(), length() Manipulating vectors and matrices reshape(), transpose operator(‘) Transpose algebra Inverse algebra Row/column vectors Matrices from vectors Semicolons – index generator Vector and matrix norms (1, 2, infinity, Frobenius) Special Matrices • Square, symmetric, skew-symmetric, diagonal, identity, upper/lower triangular, banded or Toeplitz, circulant, Vandermonde • eye(), zeros(), ones(), diag(), toeplitz(), rand(), randn(), etc. Notation and Operators • • • • Colon notation Indexing Conformance Expressions, operators (normal and elementwise), and operator precedence (+, -, *, /, ^, :, ;, ,, …, %, ‘, =, (), [], .+, .-, .*, ./, .^, .\, .’, < <=, >, >=, ==, ~=, &, |, &&, ||) • Operator Precedence Built-In Functions • cd, clc, clear, clear x, dir, pwd, exist, type, who, which • linspace • disp, length, ndims, numel, size • cross, diag, dot, end, kron, max, min, prod, reshape, sort, sum, size • det, inv, linsolve, lu, norm, null, orth, rank, rref, trace • cond – condition number and loss of precision Operators • Calculations with vectors and matrices (+, -, *, ^, .*, .^) • Operators and element-wise operators • Solving linear systems (inv, LU, QR, SVD, normal equations, pseudo-inverse) • Inconsistent • Consistent, full rank • Consistent, rank-deficient Operators • • • • Least squares fitting Robust fitting Matrix inverse Matrix functions (det, diag, eig, inv, norm, rank) Special Variables • ans, eps, i, j, NaN, pi, Inf, and –Inf Programming • Writing your own functions (implies input and output arguments, variables in local workspace) • Functions – [output_argument_list]=function function_name(input_argument_list) • Scripts (implies no input or output arguments, variables in global workspace) – Variables in the global workspace may be overwritten – Script execution can be affected by variables in the global workspace – Advisable to use functions for large and/or complicated applications Programming cont’d • • • • • • • • • Comments (%) Recursion Pretty print Operation types Sequential (commands executed in order) Conditional (if-then-else, switch-case) Iterative (for loops, do loops) Structured programing Reusable code Programming - Flow Control • Looping – for-end – Indexing (also reverse indexing) – while-end – break – continue – return Programming - Flow Control • Conditional statements – if-elseif-else-end – Relational operators (>, <, >=, <=, ==, ~=, &, |, ~) – Boolean algebra (truth tables, Karnaugh maps, DeMorgan’s Law) – bitor, bitand, bitxor, xor, and, or, not – switch-case-end Debugging • • • • • • • • • • Breakpoints Conditional breakpoints Step into Step over Return Examining variables dbquit disp() End debugging from menu Function keys Help System • Help command • Product Help Editor • Editor window • Pretty print (smart indent, ctrl-i) • Comments - % (ctrl-r, ctrl-t) I/O – Plotting and Graphics • • • • • 2D/3D Plots Plotting points and lines plot(), plot3(), surf(), mesh() Line color and style (color, linestyle) Marker color and style General I/O • • • • • To/from screen and files Binary and text file I/O format (format long, format short, etc.) Numerical format load, save, fopen, fscanf, fprintf, sprintf, fclose Image I/O • • • • imread, imwrite, imshow, imagesc, colormap Image types Image manipulation Image display (imshow, imscale) Fourier • • • • • • • Complex variables Complex algebra Complex exponential Fourier series Fourier transform Fast Fourier transform Sampling GUI Development • Event-driven programming • Design and implementation Matrix Decompositions • LU, Cholesky, QR, SVD • Condition number, accuracy, FLOPS Equation Solving • • • • • • • Matrix division inv(), \, / Pseudo-inverse LU decomposition QR decomposition SVD Operation counts versus robustness Least-Squares • • • • • • Normal equations Pseudo-inverse solution Ill-conditioning SVD solution Camera calibration 3D measurement Eigenvalues • eig() Performance • • • • Vectorization Preallocation of variables Timing (tic, toc) Profiler Toolbox Overview • • • • • Signal Processing Image Processing Statistics Optimization Parallel Processing Image Processing and Image Analysis • Morphology • Filtering • Convolution