Research: Practice on Matlab Day 5 (MIT 6.094)

Motivation is a hard thing to keep, I did finish the homework by day 5 but then I put off publishing this because my mind was preoccupied by something else. I bought a car, and now I have to fix it and practice on it so I can pass driving exam.

For now I’ll just upload without much description. I also skipped the optional homework.. Someday I’ll do it..

Resources:

Lecture 4 and Homework 4

LECTURE 4                                                                                  

Statistics

Simple histogram

close all;
%populate the scores
scores = 100*rand(1,100);

%histogram with bin centered at 5,15,...,95
hist(scores, 5:10:95);

%returns the number of occurences between specified bin edges
   %[0,10), [10,20), [30,40), ..., [90,100],
N=histc(scores, 0:10:100);

%plot barchart
bar(0:10:100, N, 'r');

Brownian motion

Simulation of Brownian random walk, then visualize it as histogram.

a=zeros(0,10000);
x=0;
for n=1:10000
    if rand()<0.5
        %move left
        x=x-1;
    else
        %move right
        x=x+1;
    end

    a(n)=x;
end
figure(1);
%plot a 50 bin histogram
hist(a,50);

Sparse matrices

Making matrices sparse so they don’t take up a lot of memory, and more efficient to process. whos details the matrices size, type and attributes.

a=zeros(100);
a(1,3)=10;
a(21,5)=pi;
whos a;
b=sparse(a);
whos b;

Data structures

Practice using cell and struct as a way to organize data

clc;
clear all;

%initializing empty cell
a=cell(1,2);
%intializing a cell
c={'hello',[1 2 3 4], rand(3,2)};
%accessing a cell element
c{1,2};

%initialize empty struct
s=struct([]);
%add fields to the struct
m(1).name = 'hell joe';
m(1).scores = [100 98 50];
m(1).year = 'G3';

people=struct('name',{'John','Mary','Leo'},...
    'age',{32,27,18},'childAge',{[2,4],1,[]});
people(3).name;

Sentence generation 

Illustrates the convenience of cell for dealing with string.

clc;
%create cell
data=cell(3,2);
data(:,1)={'andy';'warhol';'zumba'};
data(:,2)={'stinky';'smelly';'shiny'};

%display sentence
disp([data{ceil(rand*3),1} ' is ' data{ceil(rand*3),2}]);

Figure handles

Changing the figure setting using handles. Anything in the figure and plot can be customized via handles.

clc; close all;
%handle for the plotted line
L=plot(1:10,rand(1,10));
%handle for the current axis
A=gca;
%handle for current figure
F=gcf;
%current property values
get(L);
%change property values
set(A, 'FontName', 'Cambria', 'XScale', 'log');
set(L,'LineWidth',1.5,'Marker','*');

Read write images

Matlab is quite good for processing images.

%write RGB image
imwrite(rand(500,500,3),'test1.jpg');
%write indices and colormap
imwrite(ceil(rand(200)*256),jet(256),'test2.jpg');
%import image to matlab
im=imread('test1.jpg');

Making animations

The script uses movie2avi which has been deprecated, use videowriter instead.

%making animation by plotting frame by frame
    ...with pause in between
close all;
for t=1:100
    imagesc(rand(500));
    colormap(gray);
    pause(.01);
end

close all;
for n=1:100
    imagesc(rand(500));
    colormap(jet);
    M(n)=getframe;
end

%play the movie 2 times at 100fps
movie(M,2,100);

%save as avi files
movie2avi(M,'testMovie.avi','FrameRate',100,'VideoCompresssionMethod','Motion JPEG AVI');

Animated gifs

l4-4-testGif.gif

%using imwrite to save animated gif
temp=ceil(rand(300,300,1,10)*256);

imwrite(temp,jet(256),'l4-4-testGif.gif',...
    'delaytime',1/60,'loopcount',100);

Revisiting random surface as animation

This isn’t in the lecture, but I was thinking of using what I’ve learned so far by making a smooth surface using random numbers and  interp2, then have a smooth time-transition between surfaces using interp, and then customize the plot using handles such as gca, gcf, and set.

l4-revisitSurfaceGif.gif

clc;close all;

%PARAMETER
%Peaks and Valleys
N=5;
%Length
T=3;
%Graph interpolation resolution
rg=0.1;
%Time interpolation resolution
rt=0.01;
%GIF settings
gifDelay=0;
gifLoop=Inf;

%INITIALIZATION
%Time dimension for Z and sZ
t0=1:T;
ts=1:rt:T;
%meshgrid for Z0
[X0,Y0]=meshgrid(1:N,1:N);
%meshgrid for interpolated Zi
[X,Y]=meshgrid(1:rg:N,1:rg:N);
%initialize for speed
Z=zeros(length(X),length(Y),T);

%CREATE GRAPH INTERPOLATION
for i=1:T
    %Peaks and valleys
    Z0=rand(N);
    %interpolate values
    Z(:,:,i)=interp2(X0,Y0,Z0,X,Y,'cubic');
end

%SMOOTHEN TRANSITION BETWEEN GRAPH
%Permute to get interpolated dimension first:
permutedZ = permute(Z,[3 1 2]);
%create sZ for smooth transition
iPermutedZ = interp1(t0,permutedZ,ts);
%permute to return to original dimension arrangement
sZ = permute(iPermutedZ,[2 3 1]);

%ANIMATE
h=figure;
for i=1:length(ts)
    %create surface
    surf(X,Y,sZ(:,:,i))
    hold on
    contour(X,Y,sZ(:,:,i))
    hold off

    %figure settings
    colormap(jet);
    set(gca,'color','k');
    set(gcf,'color','k');
    shading interp;
    zlim([0,1]);
    %colorbar;

%             %PLOT->GIF IMPLEMENTATION 1
%             % Capture the plot as an image
%             frame = getframe(h);
%             im = frame2im(frame);
%             [imind,cm] = rgb2ind(im,256);
%             % Write to the GIF File
%             if i==1
%                 imwrite(imind,cm,'revisitSurfaceImwrite.gif','gif','DelayTime',gifDelay,'Loopcount',gifLoop);
%             else
%                 imwrite(imind,cm,'revisitSurfaceImwrite.gif','gif','DelayTime',gifDelay,'WriteMode','append');
%             end

    %PLOT->GIF IMPLEMENTATION 2, (better but still choppy)
    %add to gif
    if i==1
        gif('l4-revisitSurfaceGif.gif','DelayTime',gifDelay, ...
            'LoopCount',gifLoop,'frame',gcf);
    else
        gif
    end

    %progress
    disp([num2str(round(i/length(ts)*100),4) ' %']);
end

%view the animation
web('revisitSurfaceGif.gif');
web('revisitSurfaceImwrite.gif');

Plot ten lines

the legend function didn’t work for me for the longest time until I realized I made a variable with the same name, Matlab prioritizes variables so functions don’t work anymore.

l4-5-plotten.jpg

function l4plotTenlines
    %Exercise in debugging
    clc;
    % Click on – next to line numbers in MATLAB files
    % ¾ Each red dot that appears is a breakpoint
    % ¾ Run the program
    % ¾ The program pauses when it reaches a breakpoint
    % ¾ Use the command window to probe variables
    % ¾ Use the debugging buttons to control debugger
    clear all;
    close all;

    %make x vector
    x=0:.01:1;
    %legendNames=char(ones(10,8))

    %make a figure and plot 10 random lines, labeling each one
    figure;
    %legendNames = {};   %initialize
    for n=1:10
        plot(x,polyval(rand(3,1),x),'color',rand(1,3));
        hold on;
        legendN=['Line ' num2str(n)];
        legendNames(n,1:length(legendN))=legendN;
    end
    %label the graph
    xlabel('X');
    ylabel('Y');
    title('Ten Lines Plot');

    legendNames
    legend(legendNames);

end

HOMEWORK 4

Problem 1 – random numbers

% Homework 4 - 1
% Random variables. Make a vector of 500 random numbers from a
    % Normal distribution with mean 2 and standard deviation 5 (randn).

% After you generate the vector, verify that the sample mean and standard
    % deviation of the vector are close to 2
    % and 5 respectively (mean, std).

%I played to just see the most extreme values after 100000 runs
clc;
sigma=5;
mu=2;
v=randn(500,1).*sigma+mu;
for i=1:100000
    avg=mean(v);
    if abs(avg-mu)>abs(avgD-mu)
        avgD=avg
    end
    sd=std(v);
    if abs(sd-sigma)>abs(sdD-sigma)
        sdD=sd
    end
end

Problem 2 – flip coins

Although simple, I thought this was a good illustration to show how probability works considering the course was for lower division students.

h4-2-flipcoin.jpg

% Homework 4 - problem 2
% Flipping a coin.
% Write a script called coinTest.m to simulate sequentially flipping a coin 5000 times.
    % I will define head as 1 and tail as 0
close all;
rep=10000;
flips=round(rand(1,rep));
% Keep track of every time you get ‘heads’ and plot the running estimate of the probability of getting ‘heads’ with this coin.
headEvent=cumsum(flips);
totalEvent=1:rep;
headProbability=headEvent./totalEvent;
% Plot this running estimate along with a horizontal line at the expected value of 0.5, as below.
plot(totalEvent,headProbability,totalEvent,0.5*ones(1,rep));
% This is most easily done without a loop (useful functions: rand, round, cumsum).

Problem 3 – Poisson distribution

Plotting Poisson random numbers (poissrnd) vs the probability density function (poisspdf). Also hist is deprecated in favor of histogram.

h4-3-Poisson.jpg

% Homework 4 Problem 3
% Histogram. Generate 1000 Poisson distributed random numbers with parameter ? = 5(poissrnd).
close all;
clc;
lambda=5;
count=1000;
x=poissrnd(lambda,1,count);
xMin=min(x);
xMax=max(x);
% Get the histogram of the data and normalize the counts so that the histogram sums to 1
    ...(hist – the version that returns 2 outputs N and X, sum).
H=histogram(x,'BinEdges',[xMin:xMax],'Normalization','probability');
% Plot the normalized histogram (which is now a probability mass function) as a bar graph (bar).


% Hold on and also plot the actualPoisson probability mass function with ? = 5 as a line (poisspdf).
hold on;
index=linspace(xMin,xMax,count);
plot(index,poisspdf(lambda,index),'r');
% You can try doing this with more than 1000 samples from the Poisson distribution to get better agreement between the two.

Problem 4 – cells

Cells are interesting that they can fit anything inside one of their element. Also accessing cell element as a cell requires the use of round brackets (). While accessing cell element as their respective data types requires the use of curly brackets {}.

% Practice with cells.
% Usually, cells are most useful for storing strings, because the length of each string can be unique.
% a. Make a 3x3 cell where the first column contains the names: ‘Joe’, ’Sarah’, and ’Pat’, the
% second column contains their last names: ‘Smith’, ‘Brown’, ‘Jackson’, and the third
% column contains their salaries: $30,000, $150,000, and $120,000. Display the cell using
% disp.
employee={'Joe','Smith',30000;'Sarah','Brown',150000;'Pat','Jackson',120000};
disp(employee)

clear employee2;
employee2=cell(3,3);
employee2(:,1)={'Joe';'Sarah';'Pat'};
employee2(:,2)={'Smith';'Brown';'Jackson'};
employee2(:,3)={30000,150000,120000};
disp(employee2)

% b. Sarah gets married and changes her last name to ‘Meyers’. Make this change in the cell
% you made in a. Display the cell using disp.
employee(2,2)={'Meyers'}
% c. Pat gets promoted and gets a raise of $50,000. Change his salary by adding this amount
% to his current salary. Display the cell using disp.
employee{3,3}=employee{3,3}+50000
% The output to parts a-c should look like this:

Problem 5 – structs

Illustrates the utility of structs by using dir to list the file properties within the current work directory.

function hw4p5structs
% Using Structs. Structs are useful in many situations when dealing with diverse data. For
% example, get the contents of your current directory by typing a=dir;
A=dir;
% a. a is a struct array. What is its size? What are the names of the fields in a?
disp(['size of A: ' num2str(size(A))]);

% b. write a loop to go through all the elements of a, and if the element is not a directory,
% display the following sentence ‘File filename contains X bytes’, where filename is the
% name of the file and X is the number of bytes.

for i=1:length(A)
   disp(['File ' A(i).name ' contains ' num2str(A(i).bytes) ' bytes']);
end

% c. Write a function called displayDir.m, which will display the sizes of the files in the
% current directory when run, as below:
end

Problem 6 – handles 

Customizing plots using handles ~ gca, gcf, set, and grid

h4-5-handle.jpg

% Handles. We’ll use handles to set various properties of a figure in order to make it look like this:
% One sine wave from 0 to 2? X values (in terms of ?) Sin(X)
% a. Do all the following in a script named handlesPractice.m
% no
% b. First, make a variable x that goes from 0 to 2? , and then make y=sin(x).
res=200;
x=linspace(0,2*pi,res);
y=sin(x);
% c. Make a new figure and do plot(x,y,’r’)
figure;
plot(x,y,'r');
% d. Set the x limit to go from 0 to 2? (xlim)
xlim([0,2*pi]);
% e. Set the xtick property of the axis to be just the values [0 pi 2*pi], and set
% xticklabel to be {‘0’,’1’,’2’}. Use set and gca
set(gca,'xtick',[0,pi,2*pi],'xticklabel',{'0','1','2'});
% f. Set the ytick property of the axis to be just the values -1:.5:1. Use set and gca
set(gca,'ytick',-1:0.5:1);
% g. Turn on the grid by doing grid on.
grid on;
% h. Set the ycolor property of the axis to green, the xcolor property to cyan, and the
% color property to black (use set and gca)
set(gca,'xcolor','cyan','ycolor','black','color','black');
% i. Set the color property of the figure to a dark gray (I used [.3 .3 .3]). Use set and gcf
set(gcf,'color',[.3 .3 .3]);
% j. Add a title that says ‘One sine wave from 0 to 2?’ with fontsize 14, fontweight
% bold, and color white. Hint: to get the ? to display properly, use \pi in your string.
% Matlab uses a Tex or Latex interpreter in xlabel, ylabel, and title. You can do all this just
% by using title, no need for handles.
title('One sine wave from 0 to 2\pi', 'fontsize', 14,'fontweight', 'bold','color','white');

% k. Add the appropriate x and y labels (make sure the ? shows up that way in the x label)
% using a fontsize of 12 and color cyan for x and green for y. Use xlabel and ylabel
xlabel('multiples of \pi','fontsize',12,'color','cyan');
ylabel('sin(x\pi)','fontsize',12,'color','green');

% l. Before you copy the figure to paste it into word, look at copy options (in the figure’s Edit
% menu) and under ‘figure background color’ select ‘use figure color’.

Problem 7 – image

This question was problematic for me because at first I didn’t really understand how meshgrid really works and I flipped rows and columns (as x and y variables). and had it fine and dandy the whole time because I was dealing with square matrices/images. Until I had to deal with non-square images in this exercise. Also this exercise taught me how some image processors resizes images and results in smoothed but blurry picture. It was 2D interpolation all along.

h4-7-parthenonrgb.jpg

% Image processing
% Write a function to display a color image, as well as its red, green, and blue layers separately
function im=hw4p7image(filename)
    % The function declaration should be im=displayRGB(filename)
    % filename should be the name of the image (make the function work for jpg images only)
    clc;

    %if function is called with no argument
    if nargin == 0 %
        %filename = 'test1.jpg'      %SUCCESS: This is a square image nxn
        filename = 'https://upload.wikimedia.org/wikipedia/commons/d/da/The_Parthenon_in_Athens.jpg'  %FAIL: This is a rectangular image mxn
    end

    %this function can only accept JPG
    [pathstr,name,ext] = fileparts(filename);
    %somehow using short-circuit
    if length(ext) ~= 4 | ext ~= '.jpg'
        error('The file received is not jpg');
    end

    % im should be the final image returned as a matrix
    % To test the function, you should put a jpg file into the same directory as the function and run it with the filename
    % (include the extension, for example im=displayRGB(‘testImage jpg’))
    % You can use any picture you like, from your files or off the internet
    % Useful functions: imread, meshgrid, interp2, uint8, image, axis equal, axis tight
    imageOriginal=imread(filename);
    % a To make the program work efficiently with all image sizes, first interpolate each color layer of the original image
    % so that the larger dimension ends up with 800 pixels
    % The smaller dimension should be appropriately scaled so that the length:width ratio stays the same
    % Use interp2 with cubic interpolation to resample the image
    imageSize=size(imageOriginal);

    %convert to double
    imageOriginal=double(imageOriginal);
    %the parameter for new image
    newLength=800;

    if max(imageSize) ~= newLength
        %get all information from original image
        imageSize=size(imageOriginal);
        x=1:imageSize(2);
        y=1:imageSize(1);
        [X,Y]=meshgrid(x,y);

        %new dimension
        newSize=floor(imageSize/max(imageSize)*newLength);
        xn=linspace(1,imageSize(2),newSize(2));
        yn=linspace(1,imageSize(1),newSize(1));
        [Xn,Yn]=meshgrid(xn,yn);
        %new image
        for i=1:3
            imageResampled(:,:,i)=interp2(X,Y,imageOriginal(:,:,i),Xn,Yn,'cubic');
        end

    else
        imageResampled=imageOriginal;
    end

    %reconvert both images as uint8
    imageOriginal=uint8(imageOriginal);
    imageResampled=uint8(imageResampled);

    close all;
    %for comparison
    figure(1);image(imageOriginal);title('before resampling');
%         figure(2);image(imresize(imageOriginal,newSize));title('built in resize function');
    % b Create a composite image that is 2 times as tall as the original, and 2 times as wide
    % Place the original image in the top left, the red layer in the top right, the green layer in the bottom left, and the blue layer in the bottom right parts of this composite image
    % The function should return the composite image matrix in case you want to save it as a jpg again (before displaying or returning, convert the values to unsigned 8-bit integers using uint8)
    % Hint: To get just a single color layer, all you have to do is set the other two layers to zero
    % For example if X is an MxNx3 image, then X(:,:,2)=0; X(:,:,3)=0; will retain just the red layer
    % Include your code and the final image in your homework writeup
    imageRed=imageResampled;
    imageRed(:,:,2)=0;imageRed(:,:,3)=0;
    imageGreen=imageResampled;
    imageGreen(:,:,1)=0;imageGreen(:,:,3)=0;
    imageBlue=imageResampled;
    imageBlue(:,:,1)=0;imageBlue(:,:,2)=0;
    imageDecomposition=[imageResampled,imageRed;imageGreen,imageBlue];

    figure(3);
    image(imageDecomposition);title('image decomposition');
end

 

Problem 7 troubleshooting – interpExperiment

To troubleshoot 7. I had to realize first I might be doing something sketchy with meshgrid or interp2 so I made a smaller function to deal with those two. Voila, found my problem.

clc; clear all; close all;
%original dimensions
setOriginalSize=[3,5,3];       %param1

%random numbers
imageOriginal=rand(setOriginalSize);

%get all information from originalMat
imageSize=size(imageOriginal);
x=1:imageSize(2);
y=1:imageSize(1);
[X,Y]=meshgrid(x,y);

%new dimension
newLength=100;
newSize=floor(imageSize/max(imageSize)*newLength);
xn=linspace(1,imageSize(2),newSize(2));
yn=linspace(1,imageSize(1),newSize(1));
[Xn,Yn]=meshgrid(xn,yn);
%new image
for i=1:3
    imageResampled(:,:,i)=interp2(X,Y,imageOriginal(:,:,i),Xn,Yn,'cubic');
end
%--------------------------------------------------------------------------
figure(1);
image(imageOriginal);
figure(2);
image(imageResampled);

Problem 8 – Brownian diffusion

I thought this was pretty neat. Illustrating Brownian diffusion only using normal distribution (randn).

h4-brownianDiffusion.gif

% Animation: Brownian motion
% Write a function with the following declaration:
function hw4p8browniandiffusion(N,iter,delay,sigma)
    clc;
% The function takes in a single input N, which is an integer specifying the number of points in the simulation
    if nargin == 0
        N=100;
        iter=1000;
        delay=0.01;
        sigma=0.005;
    elseif nargin == 1
        iter=1000;
        delay=0.01;
        sigma=0.005;
    elseif nargin == 2
        delay=0.01;
        sigma=0.005;
    elseif nargin == 3
        sigma=0.005;
    end
% All the points should initially start at the origin (0,0)
    x=zeros(1,N);
    y=zeros(1,N);
% Plot all the points on a figure using ‘.’ markers, and set the axis to be
    ...square and have limits from -1 to 1 in both the x and y direction
    close all;
    P=plot(x,y,'.');
    xlim([-1,1]);ylim([-1,1]);
% To simulate Brownian motion of the points, write a 1000-iteration loop which will calculate a new x and y
    ...position for each point and will display these new positions as an animation

    for i=1:iter
% The new position of each point is obtained by adding a normally distributed random variable with a standard deviation of 0
    ... 005 to each x and y value (use randn; if you have 100 points, you need to add 100 distinct random values to
    ...the x values and 100 distinct random values to the y values)
        x=x+randn(1,N)*sigma;
        y=y+randn(1,N)*sigma;
% Each time that the new positions of all the points are calculated, plot them on the figure and pause for 0
    ... 01 seconds (it’s best to use set and the line object handle in order to
    ...update the xdata and ydata properties of the points, as mentioned in lecture)
        set(P,'XData',x,'YData',y,'Marker','.');
        pause(delay);
    end
% What you will see is a simulation of diffusion, wherein the particles randomly move away from the center of the figure
% For example, 100 points look like the figures below at the start, middle, and end of the simulation:
end

Optional Problem – Julia Set animation

 

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google+ photo

You are commenting using your Google+ account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s