ESS116
ESS 116
Introduction to Data Analysis in Earth Science
Image Credit: NASA
Instructor: Mathieu Morlighem
E-mail: mmorligh@uci.edu (include ESS116 in subject line)
Office Hours: 3218 Croul Hall, Friday 2:00 pm – 3:00 pm
This content is protected and may not be shared uploaded or distributed
Lecture 7 quick review
Lecture 8 – Image processing
Images in MATLAB
Images indexed to a color map
True color images (RGB)
Examples of Image processing
Midterm exam solutions
Today’s lecture
Review Lecture 7
Polynomials in MATLAB
Any polynomial can be represented as a row vector:
The length of the row vector indicates the degree of the polynomial
degree of the polynomial = length(coeff) – 1
Use polyval to compute the value of a polynomial for x
Review Lecture 7
>> coeff = [3 2.7 1 -5.7];
>> polyval(coeff,2.5)
ans =
60.5500
>> coeff = [3 2.7 1 -5.7];
>> x=0:0.2:0.6;
>> polyval(coeff,x)
ans =
-5.70 -5.368 -4.6760 -3.48
Which of the following represents a linear polynomial?
coeff = 0;
coeff = 3;
coeff = [0 3];
coeff = [3 0];
both C and D
i>Clicker question
Fitting Data with a polynomial
Use polyfit
For nonlinear data, do not go higher than 4 or 5 (do not over-fit)
Interpolation:
Review Lecture 7
>> coeff = polyfit(datax,datay,1);
Use polyfit or interp1 (linear or spline)
Never extrapolate…
>> xI=1:5:100; %We want to interpolate on these points
>> coeff = polyfit(datax,datay,1);
>> yI1=polyval(coeff ,xI); %Polynomial interp
>> yI2=interp1(datax,datay,xI,’linear’); %Linear interp
>> yI3=interp1(datax,datay,xI,’spline’); %Spline interp
Lecture 8 – Image Processing
Earth science is a very visual discipline
Graphs, Maps, Field Photos, Satellite images
Flood detection, ice sheet retreat, forest fires,…
Because of this, all Earth scientists should have:
Basic knowledge about graphics file types
Basic idea on how to process images
Image Basics
Photomicrograph of peridotite (mantle rock)
NED Dataset
Sediment Core
Landsat image of Death Valley, CA
Cathodoluminescence image of granite
Raster Image:
Made of pixels (each pixel has a uniform color)
When you scale it, quality changes (finite number of pixels)
Common formats: jpg, png, gif, tif, bmp
Vector Image:
Made of vector objects (maths, not pixel-based)
Can scale to any size
Common Formats: ai, eps, ps, svg, wmf, pdf
Many Formats can contain both vector and raster images
For example, ai, pdf, eps, ps, wmf
How can you tell? Zoom in on the image
Which is better?
For graphs/plots: typically vector (unless dataset is HUGE)
For photos: raster
Image File Types: Raster vs Vector
Raster vs Vector
Raster (jpeg)
Vector (pdf)
Image Processing in MATLAB
MATLAB provides functions that read raster images
Each pixel’s value (color) is stored as a number in a matrix
MATLAB also has an image processing toolbox with TONS of image processing options
In this class, we will only use the image functions that are part of the standard MATLAB libraries
Image Processing in MATLAB
Novarupta Volcano, Aleutian Islands, Alaska
Purples/Reds: Volcanic ash from 1912 Eruption
Blues: Snow/Ice
Mississippi River meanders & oxbows near Memphis, TN
Yukon River delta, Alaska
Landsat false-color image examples…
Before we start processing images, we need to talk about how computers represent colors as numbers
Three common color models
RGB: An additive model
works like light
CMYK: A subtractive model
works like ink
HSV: A cone-shaped model
useful for shading colors
Color Models
For variables of class uint8 (and 8-bit image)
0-255 are the possible integer values (same as the uint8 class!)
0 is minimum for any RGB color
255 is max for any RGB color
One color: must specify Red (R), Green (G) and Blue (B) values
Grayscale images only need one value (0=black, 255=white)
The RGB Color Model
[ 0 0 0 ]
[ 255 255 255 ]
[ 75 75 75 ]
[ 200 200 200 ]
[ 255 0 0 ]
[ 100 0 0 ]
[ 0 255 0 ]
[ 0 100 0 ]
[ 0 0 255 ]
[ 0 0 100 ]
[ 255 255 0 ]
[ 0 255 255 ]
[ 255 0 255 ]
[ 237 125 49 ]
Total number of different colors: 2563 = 16,777,216
Just to annoy us, MATLAB requires colormap RGB values to be values between 0 and 1 for variables of class double
0 is minimum for any RGB color
1 is max for any RGB color
One color: must specify Red (R), Green (G) and Blue (B) values
To convert from the 0-255 system, cast to double and divide by 255
Grayscale images really only need one value (0=black, 1=white)
Colormaps in MATLAB
[ 0 0 0 ]
[ 1 1 1 ]
[ 0.29 0.29 0.29 ]
[ 0.78 0.78 0.78 ]
[ 1 0 0 ]
[ 0.39 0 0 ]
[ 0 1 0 ]
[ 0 0.39 0 ]
[ 0 0 1 ]
[ 0 0 0.39 ]
[ 1 1 0 ]
[ 0 1 1 ]
[ 1 0 1 ]
[ 0.93 0.49 0.19 ]
MATLAB can represent color in images in two basic ways
True Color or RGB
The three color components are stored in a m x n x 3 matrix. i.e. a 3D matrix.
(:, :, 1) R-values; (:, :, 2) G-values; (:, :, 3) B-values
Indexed to a Colormap
Colors are stored as a single integer value that corresponds to a row in a colormap matrix.
Colormap stores the RGB values
MATLAB & Images
Image Matrix
Colormap Matrix
4 Pixel Image
Images indexed to a colormap
‘image’ plots a matrix as an image
Don’t forget to specify the colormap (If not, you get the default 64 color parula colormap)
Images Represented as Colormaps
MATLAB provides several built-in colormaps
The command colormap is rather useful
It can set the current colormap
Can be a built-in map or a custom n x 3 matrix
Built-in colormaps can be easily scaled
E.g. “jet(256)” returns a 256 x 3 matrix that follows the color scheme of the built-in “jet” colormap.
Built-In Colormaps
We would like to create a French flag using
iMat = [1 2 3]
What colormap should be used?
cmap = [1 0 0;0 1 0;0 0 1];
cmap = [1 0 0;1 1 1;0 0 1];
cmap = [0 0 1;1 1 1;1 0 0];
cmap = [1 0 0;0 0 0;0 0 1];
Don’t know…
i>Clicker question
If you exceed the color map values, you get either the min or max color
Exceed Max/Min Indexed Color?
imagesc will scale the image matrix to use the full range of colors
Exceed Max/Min Indexed Color?
“True color” Images (RGB)
If ‘image’ is passed a 3D matrix
It is assumed to be a true color image
1st level R-values
2nd level G-values
3rd level B-values
If ‘image’ is passed a 2D matrix
It is assumed to be a colormap indexed image
True Color Images
If we would like to create a French flag using an RGB image, which of the following 3D matrices would you use?
iMat(:,:,1) = [0 1 1];
iMat(:,:,2) = [0 1 0];
iMat(:,:,3) = [1 1 0];
iMat(:,:,1) = [0 0 1];
iMat(:,:,2) = [1 1 1];
iMat(:,:,3) = [1 0 0];
Don’t know…
i>Clicker question
iMat =
i>Clicker question
[ 1 1 0 ]
[ 0 1 0 ]
[ 0 1 1 ]
Blue for all pixels
Green for all pixels
Red for all pixels
BLUE
WHITE
RED
Why does MATLAB offer two ways to store image colors?
Flexibility. It is always good to give users options
Grayscale vs Color Images
These images are typically read in MATLAB in different ways
Grayscale: Indexed Color
Color image: True Color or RGB image
Why have both Types of Images?
Dendrites example (indexed)
Reading in Grayscale Images
‘imread’ can read in most standard grayscale raster image types
.jpg .gif .png, etc…
Stores image as a rectangular matrix
Each entry represents one pixel’s grayscale value
0-255 (black to white)
Let’s read this image into MATLAB
The image is 814 x 531 pixels
Matrix will be 531 x 814
Try to automate detection of the dendrites
Reading in Grayscale Images
Cropped grayscale image of dendrites
Increased contrast (using Photoshop)
If ‘image’ is given a 2D matrix, it assumes the image is indexed color
Reading in Grayscale Images
minGray=75;
imageMat = imread(‘den3_highContrast.jpg’);
%Pre-Allocate to ones (white pixels)
denMask = ones(size(imageMat));
%Loop through the rows and cols of the image
[rows cols]=size(imageMat);
for i=1:rows
for j=1:cols
if imageMat(i,j)
waterMask(i,j,1)=0;
waterMask (i,j,2)=0;
waterMask (i,j,3)=0.70;
end
end
end
%plot result
subplot(2,1,1);
imagesc(iMat);
axis(‘equal’,’tight’);
colormap(gray);
title(‘RAW image’)
subplot(2,1,2);
imagesc(waterMask);
axis(‘equal’,’tight’);
title([‘Identified Water pixels, Blue Value > ‘ num2str(minBlue)]);
Reading in True Color / RGB Images
Are all of the selected pixels actually blue? (No)
How could I improve water detection?
How could I calculate water surface area?
MATLAB Commands to remember
Lab 8: Image processing
Lecture 9: time Series
What’s next?
3×3 + 2.7×2 + x� 5.7
⇥
3 2.7 1 �5.7
⇤
UCI14_EarthSysSci_2Lines_Blue
2
66
4
1 0 0
0 1 0
0 0 1
1 1 1
3
77
5
1 3
2 4
�
ESS116: MATLAB Cheat Sheet
1 Path and file operations
cd Change Directory (followed by absolute or relative path of a directory)
cd ../../Shared (relative path)
cd /Users/Shared (absolute path)
pwd display current directory’s absolute path (Path Working Directory)
ls display list of files and directories in the current directory
(can be followed by a path and/or file name pattern with *)
ls ../file*mat
ls *.txt
ls /Users/mmorligh/Desktop/
copyfile copy existing file into a new directory, and/or rename a file
copyfile(‘/Users/Shared/foo.txt’,’.’);
copyfile(‘foo.txt’,’bar.txt’);
mkdir create a directory
mkdir Lab1
2 Fundamental MATLAB classes
double floating point number (1.52, pi, …) → MATLAB’s default type
int8 Integer between -128 and 127 (8 bits, saves memory)
uint8 Unsigned integer between 0 and 255 (used primarily for images)
int16 Integer between -32768 and 32767 (16 bits)
logical true/false
string data type for text (str = ‘This is a string’;)
cell cell array, used by textscan
3 Matrices
Use square [] to create a matrix, and ; to separate rows
A=[1 2 3;4 5 6;7 8 9];
ones, zeros create a matrix full of ones or zeros
A=ones(5,2);
‘ transpose a matrix
B=A’;
length return length of a vector (do not use for matrices)
size returns the size of a matrix
(number of rows then columns, then 3rd dimension if 3D, etc)
[nrows,ncols]=size(A); [nrows,ncols,nlayers]=size(A3D);
linspace and : to create vectors
A=2:3:100;
A=linspace(2,100,10);
find return the linear indices where a condition on the elements of a matrix is met
pos=find(A==−9999);
pos=find(A>100);
Extract the first 10 even columns of a matrix
B=A(:,2:2:20);
Removing elements: use empty brackets
A(:,2)= [];
Concatenate matrices
A=’This is ‘; B=[A ‘an example’];
Replacing elements in a matrix (use either linear or row,col notation)
A(10,3)=5.5;
pos=find(A==−9999);
A(pos)= 0;
Element-by-element operation: use a dot (.) before the operator
A= C.*D;
4 I/O
load loads a MATLAB file (*.mat) into the workspace, or a text file with only numbers
and consistent number of columns
load(‘data.mat’);
data=load(‘data.txt’);
textscan loads a text file into a cell array (as many elements as there are columns in the file)
Use %d for integers, %f for floating point numbers %s for strings
fid = fopen(‘filename’);
data = textscan(fid,’%d %f %s %s’,’Headerlines’,5);
fclose(fid);
%Put first column in A, and second column in B
A = data{1}; B = data{2};
5 fprintf
fprintf print text (and variables) to the screen. First argument is a string with placeholders.
fprintf(‘The radius is %7.2f and A = %d !!\n’,EarthRadius,10);
– Special characters: \n (new line) %% (percent sign) ” (apostrophe)
– Variable specifiers: %s (string) %d (integer) %e (exponential) %f (float)
– %010.3f: leading 0, 10 total spaces, 3 decimals. Ex: 000003.142
6 Visualization
plot displays a list of points (x,y)
plot(x,y,’−r’);
plot(x,y,’r+:’,’MarkerFaceColor’,’g’,’MarkerSize’,5,’LineWidth’,2);
axis controls x and y axes
axis([xmin xmax ymin ymax]);
axis equal tight
legend adds a legend to previously plotted curves
legend(‘First curve’,’second curve’)
figure creates a new figure window
figure(2)
xlabel/ylabel/title control x/y axis labels and plot title
xlabel(‘Distance (km)’);
hold on keep current plot so that whatever follows is plotted on the same plot
subplot divide figure into several subplots
subplot(2,3,1)
histogram make a histogram for a vector
histogram(tmax,20);
histogram(tmax,round(sqrt(length(tmax))));
7 Relational and Logical operators
== equal to, ~= not equal to, > greater than, >= greater than or equal to,
< less than, <= less than or equal to.
&& and, || or, ~ not.
A=( (1>10)|| (3~=4));
8 If/elseif/else and for loops (examples)
Counting algorithm
%Initialize counters
counter1 = 0;
counter2 = 0;
%Go over all of the elements of T and increment counters when a condition is met
for i=1:length(T)
if T(i)>100
counter1 = counter1+1
elseif T(i)<0
counter2 = counter2+1
end
end
fprintf('Found %d days with T>%f, and %d with T<%f\n',counter1,100,counter2,0)
Extracting (after counting!)
%You first need to count how many times T>100 (for example), then: allocate memory
hotdays = zeros(counter1,1);
%Go through T, again, and store temperatures>100 in hotdays
count = 1;
for i=1:length(T)
if T(i)>100
hotdays(count)=T(i);
count = count+1;
end
end
9 Statistic
mean computes mean of a vector
median computes median of a vector
std computes standard deviation
min returns minimum value in a vector
max returns minimum value in a vector
skewness returns skewness
kurtosis returns kurtosis
normcdf/tcdf/chi2cdf cumulative density function for a normal, t and χ2 distributions
norminv/tinv/chi2inv inverse of the cumulative density function
p0 = normcdf(x0,mu,sigma);
x0 = norminv(p0,mu,sigma);
p0 = tcdf(x0,V);
x0 = tinv(p0,V);
10 Polynomials and interpolations
polyval returns the value of a polynomials (represented by its coefficient) for some x
coeff = [3 2.7 1 −5.7];
x=0:0.2:0.6;
y=polyval(coeff,x)
polyfit returns the coefficient of the polynomials that best fit data points
coeff = polyfit(datax,datay,3); %3 means cubic polynomial
interp1 interpolates between data points (spline or linear)
y1=interp1(datax,datay,’linear’);
y2=interp1(datax,datay,’spline’);
11 Image processing
imread loads an image (as a matrix) into the workspace
A=imread(‘image.png’)
image display an image (2D or 3D)
imagesc display a 2D image and scale indices to use all the colors in the color map
colormap set a colormap (only for indexed images)
Indexed Images (2D)
Indexed image, need two matrices:
iMat a 2D matrix with indices
cMap is a nx3 2D matrices with the RGB code for each index
To display this image:
image(iMat);
colormap(cMap);
If the colormap is not consistent with indices, you need to use imagesc(iMat).
True-color image (3D, RGB)
No need to prescribe a colormap:
iMat(:,:,1) Red matrix (between 0–255 if uint8, or 0–1 if double)
iMat(:,:,2) Green matrix
iMat(:,:,3) Blue matrix
To display this image: image(iMat).
12 Miscellaneous
rand returns a random floating point number between 0 and 1
x=rand
x=rand(10,2)
round round input to closest integer
A=round(rand*10);
whos displays list of all variables in MATLAB workspace
tic/toc displays cpu time for a chunk of code
sqrt square root
13 Functions
Calling a function: [output1,output2] = functionname(arg1,arg2);
Function header (top lines of the file that implements this function):
function [output1,output2] = functionname(arg1,arg2)
% H1 line: describe what the function does
ESS116, M. Morlighem, Updated: March 14, 2019