MATLAB is a powerful tool, and not difficult to learn by yourself. Instead of let me explain all the details, you should learn the basics on your own. Please go through the assigned readings before starting this lesson. Please remember, you will forget everything soon if you do not practice at all.

§ Reading

Official MATLAB guide: chapters 1 to 4, 7
Interactive MATLAB Tutorials: Tutorials 1 to 8

§ Covered contents

File format: regular .m / function .m / .fig / .mat
Path
Variables:

local variable / global variable
Multidimensional arrays:

Multidimensional Matrices: M = rand(4,4,3);

Multidimensional Cell Arrays: S{1} = ‘hello’; S{2} = rand(4,4,3);
Multidimensional Structure Arrays: A.b = rand(4,4,3);

Flow control:

if, else

if a < 30

disp('small')

elseif a < 80

disp('medium')

else

disp('large')

end

switch

switch dayString

case 'Monday'

disp('Day 1')

case 'Wednesday'

disp('Day 3')

otherwise

disp('Other')

end

for loop
while loop: try not to use while loop!
continue: passes control to the next iteration of the for loop or while loop in which it appears, skipping any remaining statements in the body of the loop
break: terminates the execution of a for or while loop. Statements in the loop that appear after the break statement are not executed
return: terminates the current sequence of commands and returns control to the invoking function or to the keyboard.

Built-in functions

help / lookfor / doc
clear / clc / clf / control+c

Simple 2D plots

plot(x,y,’r-*’)
title / xlabel / ylabel / xlim / ylim
legend
gcf / gca
line / rectangle
colors: (w)hite, (y)ellow, (m)agenta, (r)ed, (c)yan, (g)reen, (b)lue, or blac(k)
line styles: -, --, :, -.
markers: +, o, *, ., ^, v, >, <, x, s, d, p, h

Exporting plots

vector format (pdf) / bitmap format (bmp, png, tif)

§ Examples

(1) Multiplication table

clear all; close all; clc

for ii=1:9

for jj=1:9

table1(ii,jj) = ii*jj;

end

table2 = [1:9]'*[1:9];

(2) 2D Plot

clear all; close all; clc

count = 0;

for ii=0:0.1:10

count = count + 1;

x1(count) = ii;

y1(count) = sin(ii);

end

x2 = 0:0.1:10;

y2 = cos(x2);

x3 = linspace(0,10,101);

y3 = sin(x3).*cos(x3);

figure;

set(gcf,'position',[50 50 600 400])

hold on

h1 = plot(x1,y1,'b-');

get(h1)

set(h1,'linewidth',2)

plot(x2,y2,'ro--')

plot(x3,y3,'m*')

hold off

title('Wave')

xlabel('time (\mu s)')

ylabel('displacement (mm)')

grid on

box on

xlim([0 11])

ylim([-1 1]*1.5)

h = legend({'y=sin(x)';'y=cos(x)';'y=sin(x)*cos(x)'});

set(h,'location','NorthEast')

get(gca)

set(gca,'xtick',[0:2.5:10])

line([0 10],[-1.25 -1.25],'Linewidth',2,'Color',[0 0.5 0])

text(8,-1.35,'This is a line')

h1 = rectangle('position',[2.5 1 2.5 0.5],'Curvature', [1 1]);

set(h1,'FaceColor','w','EdgeColor','r')

text(2.75,1.25,'This is an ellipse')

§ Exercise

Ray tracing simulation

Given:

Snell’s law: n1*sin(theta1) = n2*sin(theta2)

n_air = 1;

n_water = 1.33;

n_PDMS = 1.41;

fiber OD: 155 um

fiber ID: 50 um

emitting angle: 15 degrees w.r.t. horizontal

line equation: y = ax + b

ellipse equation: (x-h)^2/a^2 + (y-k)^2/b^2 = 1

dot product: a dot b = a*b*cos(theta)

Goal: design a collimating lens as shown in figure below

Lesson 2: Advanced graphs

Please remember, you will only learn from practice.

§ Reading

Official MATLAB guide: chapters 1 to 4, 7
Interactive MATLAB Tutorials: Tutorials 3 to 9

§ Covered Contents

Advanced graphs

subplot: plot multiple graphs in one figure window
2D graphs: plot / bar / stairs / errorbar / stem / polar → example 1
3D graphs: mesh / surf / surfl / contour / contour3 → example 2
images: imagesc → example 3
graph controls

colormap
caxis
colorbar

Saving, Importing, and Extracting Data

Importing data/image

use Import Wizard
use importdata function to load data: A = importdata(‘filename.txt’);
use imread function to read image: A = imread(‘filename’,’extension’);

Extracting data: row, column indexing / B = A(2:4,3:4);
Saving data/figure:

use save to save workspace variables to disk
use saveas or print to save a figure in desired output format and resolution

Basic curve fitting

Curve fitting with polynomials

use built-in Basic Fitting GUI
use polyfit and polyval functions → p = polyfit(x,y,n) finds the coefficients of a polynomial p(x) of degree n that fits the data

Ray tracing

exercise from Lesson 1
spherical aberration
weighting + multiple source

§ Examples

(1) 2D plots

function H = Plots_2D(x,y)

clc, clf

% "nargin" returns the number of input arguments that were used to call a

% user-defined function. If zero, the default parameters will be used.

if nargin==0

x = 0:0.1:4;

y = sin(x.^2).*exp(-x);

e = rand(size(x))/40;

end

H.a1 = subplot(2,3,1);

H.h1 = plot(H.a1,x,y,'r.-','linewidth',2);

Sub_axis_setting(H.a1,'plot','x','y',[0 4],[-0.15 0.35])

H.a2 = subplot(2,3,2);

H.h2 = bar(H.a2,x,y);

Sub_axis_setting(H.a2,'bar','x','y',[0 4],[-0.15 0.35])

H.a3 = subplot(2,3,3);

H.h3 = stairs(H.a3,x,y);

Sub_axis_setting(H.a3,'stairs','x','y',[0 4],[-0.15 0.35])

H.a4 = subplot(2,3,4);

H.h4 = errorbar(H.a4,x,y,e);

Sub_axis_setting(H.a4,'errorbar','x','y',[0 4],[-0.15 0.35])

H.a5 = subplot(2,3,5);

H.h5 = stem(H.a5,x,y);

Sub_axis_setting(H.a5,'stem','x','y',[0 4],[-0.15 0.35])

H.a6 = subplot(2,3,6);

H.h6 = polar(H.a6,x,y);

Sub_axis_setting(H.a6,'polar')

return

function Sub_axis_setting(a,mytitle,myxlabel,myylabel,myxlim,myylim)

% Subfunciton inside a user-defined function can be called by the user-defined function only.

myrun = {

'set(a,''fontsize'',12)'

'title(a,mytitle,''fontsize'',12)'

'xlabel(a,myxlabel,''fontsize'',12)'

'ylabel(a,myylabel,''fontsize'',12)'

'xlim(a,myxlim)'

'ylim(a,myylim)'

};

if nargin==0

disp('No input arguments!')

return

end

for ii=1:nargin

eval(myrun{ii})

end

return

(2) 3D plots

clear all; clc, clf

[X,Y] = meshgrid([-2:.25:2]);

Z = X.*exp(-X.^2-Y.^2);

% mesh plot

subplot(2,3,1);

mesh(Z);

colormap(hsv)

title('mesh','fontsize',12)

% surface plot

subplot(2,3,2);

surf(Z);

colormap(jet)

title('surf','fontsize',12)

% surface plot with light and shading

subplot(2,3,3);

surfl(Z);

shading interp;

colormap(summer)

title('surfl','fontsize',12)

% contour plot

subplot(2,3,4);

contour(Z);

colormap(hsv)

title('contour','fontsize',12)

% 3D contour plot

subplot(2,3,5);

contour3(Z);

colormap(hsv)

title('contour3','fontsize',12)

% 3D contour + surface plot

subplot(2,3,6);

contour3(X,Y,Z,30)

surface(X,Y,Z,'EdgeColor',[.8 .8 .8],'FaceColor','none')

grid off

view(-15,25)

colormap(cool)

title('contour3 + surf','fontsize',12)

return

(3) Plot images

clear all; close all; clc

figure('position',[50 500 800 400])

% load and plot image

a1 = subplot(1,3,[1 2]);

Z = imread('example3_input','png');

%Z = importdata('example3_input.png');

imagesc(Z)

axis image

colormap(gray)

title(['image resolution: ' num2str(size(Z,1)) ' x ' num2str(size(Z,2))])

xlabel('x position (pixel)')

ylabel('y position (pixel)')

% graph controls

colormap(jet)

colorbar

zmax = 200;

caxis([0 zmax])

% draw a cut line

cutx = 237;

line([1 1]*cutx,[1 size(Z,1)],'Color','w','LineStyle','--')

% plot the light intensity along the cut line

a2 = subplot(1,3,3);

x = 1:size(Z,1);

y = Z(:,cutx);

plot(a2,x,y)

title(['light intensity along y=' num2str(cutx)])

xlabel('x position (pixel)')

xlim([x(1) x(end)])

ylim([0 zmax])

% save figure

P1 = get(gcf,'PaperPosition');

P2 = get(gcf,'Position');

set(gcf,'PaperPosition',[P1(1:3) P1(3)*(P2(4)/P2(3))]);

saveas(gcf,'example3_imagesc_saveas','png')

print -dpng -r300 'example3_imagesc_print' % save as png file with 300 dpi resolution

print -dpdf -r300 'example3_imagesc_pdf' % save as pdf file with 300 dpi resolution

return

Lesson 3: Useful functions

§ Covered contents

Save data to file

use save to save variables to tab-separated text file → example 1

Import data from a file and grid data

textread: read data from text file; write to multiple outputs → example 2
griddata: fits a surface of the form z = f(x,y) to the data in the (usually) nonuniformly spaced vectors (x,y,z) → example 2

Curve fitting with user-defined error function

Gaussian function: y = a*exp(-(x-b)2 / (2*c2)) → example 3
Parabolic function: np2 = n02(1-alpha2 x2) → example 3
Hyperbolic secant function: nh2 = ns2 + (n02-ns2) sech2(alpha*x) → example 3
curve fitting: use fminsearch function (Find minimum of unconstrained multivariable function using derivative-free method) → example 4

Band pass filtering

read my blog: http://scslin.blogspot.com/2011/06/matlab-low-pass-filtering.html

§ Examples

(1) Save variables to file

clear all; close all; clc

% Sample a function at 100 random points between ±2.0

x = rand(100,1)*4-2;

y = rand(100,1)*4-2;

z = x.*exp(-x.^2-y.^2);

plot3(x,y,z,'.')

grid on

data = [x y z];

% save sample to a text file in 8-digit ASCII form

save 'sample.txt' data -ASCII

(2) Import and grid data

clear all; close all; clc

% import data from a file

[x y z] = textread(...

'sample.txt',... % Location of the file

'%f %f %f',... % Format: read a floating-point value

'delimiter',' ',... % Act as delimiters between elements

'headerlines',0 ... % Ignores lines at the beginning

);

% plot data

plot3(x,y,z,'o')

grid on

data = [x y z];

% grid data

[X Y] = meshgrid(-2:0.1:2);

Z = griddata(x,y,z,X,Y);

hold on

mesh(X,Y,Z)

(3) Gaussian, parabolic, and hyperbolic secant functions

clear all; close all; clc

% Gaussian function

a = 1.0; % a is the height of the Gaussian peak

b = 0.2; % b is the position of the center of the peak

c = 0.5; % c controls the width of the "bump"

x = linspace(-2,2,101);

y = a*exp(-(x-b).^2 / (2*c^2));

subplot(1,3,1)

plot(x,y,'r.-')

title('Gaussian function')

% Parabolic function

n0 = 1.445; % n0 is the height of the Parabolic peak

alpha = 0.2; % alpha is the gradient coefficient

np = n0.*sqrt(1-alpha^2.*x.^2);

subplot(1,3,2)

plot(x,np,'b.-')

title('Parabolic function')

% Hyperbolic secant function

n0 = 1.445; % n0 is the height of the Parabolic peak

ns = 1.2; % ns is the cutoff value

alpha = 0.2; % alpha is the gradient coefficient

ns = 1.2;

nh = sqrt(ns^2 + (n0^2-ns^2).*sech(alpha.*x).^2);

subplot(1,3,3)

plot(x,nh,'g.-')

title('Hyperbolic secant function')

(4) Curve fitting with user-defined error function

clear all; close all; clc

% Parabolic function

n0 = 1.445; % n0 is the height of the Parabolic peak

alpha = 0.2; % alpha is the gradient coefficient

x = linspace(-2,2,51);

np = n0.*sqrt(1-alpha^2.*x.^2);

plot(x,np,'b.')

% fitting with a Gaussian function

guess = [1.4 0 0.5]; % initial guesses: a b c of a Gaussian function

Vars = fminsearch('fun_el',guess,[],[x' np']);

ng = Vars(1)*exp(-(x-Vars(2)).^2 / (2*Vars(3)^2));

hold on

plot(x,ng,'r')

legend('Parabolic function','Gaussian fit')

Below is a separate .m file

function E = fun_el(g, data)

x = data(:,1);

y = data(:,2);

% fit with a Gaussian function

model_y = g(1)*exp(-(x-g(2)).^2 / (2*g(3)^2));

% error function

E = sum((y-model_y).^2);

return