OTGA UEM RTC Introduo ao Matlab introduo de

OTGA UEM RTC Introdução ao Matlab: introdução de dados, batimetria, linhas de contorno, traçar diagramas e perfis 24 Janeiro 2017 Fialho Nehama, UEM-ESCMC, Quelimane Moçambique With the support of the Government of Flanders, Belgium

Matlab – short for matrix library (matrices are the building blocks of Matlab. Matrices (or arrays) can be of all dimensions: Simple elements – 1 x 1 matrices, vectors – 1 x. N matrices, mulitple dimensional matrices (Mx. Nx…). Matlab has many built in functions and hundreds of other functions in ‘toolboxes’. For example – the statistics toolbox, containing functions for calculating correlations, mean values, variances……. . There are many matlab tutorials you can download on the internet, This is a combination of ones written by Tore Furevik and Seb Swart and Bjorn Backeberg

Change visible details Current directory Command window

Arrows on key board are useful shortcuts in Matlab. They are useful for flicking through your recent commands. The up arrow displays the last command – sin(1). Pressing it twice will bring up the ‘why’ command. The command can be edited by pressing the left arrow once and deleting 1 and putting 2. Can also put the first few letters of a command press the tab button Or a few letters of a previous command the up arrow will fill it in. Matlab is CASE SENSITIVE…. try SIN(2) If you get into a struggle – ‘ctrl c’

Matlab files are saved as. m files

Matlab file formats: . m (programs). mat (saved data). fig (saved figures Can specify to save figs as anything else, eg. jpeg, . tiff…… Can import and export data in any file (eg netcdf)

At the >> prompt, one can directly enter commands. As an example, this is one way to compute 2+2: >> a = 2 a=2 >> b = 2 b=2 >> c = a + b c=4 It is often annoying to have Matlab always print out the value of each variable. To avoid this, put a semicolon after the commands: >> a = 2; >> b = 2; >> c = a + b; >> c c=4 Only the final line produces output. Semicolons can also be used to string together more than one command on the same line: >> a = 2; b = 2; c = a + b; c c = 4

Of course Matlab will also allow more complicated operations: >> a = 2; >> b = -12; >> c = 16; >> quad 1 = (-b + sqrt(b^2 - 4*a*c)) / (2*a) quad 1 = 4 quad 2 = (-b - sqrt(b^2 - 4*a*c)) / (2*a)

The lookfor command is also very useful. It searches through matlab functions to look for whether matlab has a function for the desired purpose. Try >> lookfor mean It only looks through the first line, to search the entire function: >>lookfor mean –all Or if you have an idea already… >>which mean

>>gheko=[1, 2, 3; 4, 5, 6; 7, 8, 9] >>whos dog

Make the following matrix: X=[1 0 4 9 2; 3 9 4 2 5; 3 4 2 6 4; 2 3 8 4 3; 2 9 3 2 11]; What is the mean value of each column and row? What is the minimum value of each column What is the mean, minimum and standard value and variance of the entire matrix Type mean=2 and try and calculate the mean again? How can you get rid of this? In Matlab there are several predefined constants (as there are functions). These include pi, i and j. The first being the ratio between the circumference and diameter of a circle and the latter 2 are both the imaginary unit vector What is the area and volume of a sphere with radius 6371 km What is exp(i*pi)? Calculate (1+i)*4 and sin(2)+i*sinh(2 i)

>> mean(X, 1) ans = 2. 2000 5. 0000 4. 2000 4. 6000 5. 0000 >> mean(X, 2) ans = 3. 2000 4. 6000 3. 8000 4. 0000 5. 4000 >> min(X, [], 2) >>mean(X(: )) • Circle Area = π • r² = ¼ • π • d² • Sphere Volume = 4/3 • π • r³ = ( π • d³)/6 • Euler's formula: e^(i pi) = -1

What about nansum and nanmean

Creating matrices: >>A=3. 14 >>B=[1 9 3 5 2 4 6 8 3 9] 1 x 10 To make a more general matrix, separate by; >>C=[1 2 3 4; 3 4 5 6; 8 6 4 2] To extract elements: C(1, 1) C(2, 3) Matlab interprets a single index number as a row number. If the index number is larger than the number of rows, it will start over in the next column. C(5) becomes 4 and C(8) becomes 5. Rand(2, 3) Ones(2, 5) Magic(5) – what is special about this?

Colons!!!! >a=1: 10 >>a=1: 2: 10; >>b=9: -3: -9; >>c=0: 0. 2: 2; >>d=0: pi/6: 2*pi; >>C=[1: 4; 3: 6; 8: -2: 2]; >>C(1: 2, 3) >>C(3, 2: 4) >>C(: , 3: 4) Change sequence: >>D=C(: , [2 1 3 4]) interchanges the first and second column. Or to remove 3 rd column >>C 1=C(: , [1: 2 4]);

Or join matrices: >>C 2=[C C C 1; C 1 C C]; Or else C(3: end); get all elements from 3 rd to last (if we don’t know length) Matrix dimensions >>x=rand(2, 3, 7, 5); >>size(x) >>prod(size(x)); >>length(x) >>max(size(x)); Change the matrix dimensions without altering the sequence: >>y=reshape(x, 2, 3, 35); size(y); >>z=reshape(x, prod(size(x)), 1); size(z)

Reshape doesn’t change elements. Useful if eg have a matrix with values for every month over a period of years and we want the mean value for just January >>[val, ind]=max(x(: )) %find the max value of the matrix and it’s position in the sequence (between 1 and 210). Although it is usually more helpful to know the position in the matrix >>[a, b, c, d]=ind 2 sub(ind, size(x)) %find the max value in x(a, b, c, d) we need 4 indices since we have a 4 D matrix Check with help function >>[a, b, c, d]=ind 2 sub(size(x), ind) Testran program!!!!! Linefit program!!!!

>>x=[1: 5: 30 Na. N 20: -2: 4 Inf]; >>I=isfinite(x) %returns a vector where the elements of x are finite numbers and 0 where they are not…. find the indices to all elements of value 1 for I >>J=find(I) %elements 7 and 17 are omitted to remove them >>x=x(J) Alternately >>x=x(isfinite(x)) Replace all numbers >15 with 0 >> i=find(x>15); >> x(i)=0 >>x(find(x>15))=0; OR >>x(isprime(x))=x(isprime(x)). ^2 %squares all prime numbers >>y=x(find(x>=10) %create a new vector with all numbers >10 isnumeric and ischar checks if a variable is a number or character A=‘word’ and B=100; isnumeric(A)………

Generate a sequence of numbers from 0 to 10 with step length 0. 01 Calculate the sin and cosine to all the numbers and visualize with plot >>x=0: 0. 01: 10; >>y 1=sin(x); >>y 2=cos(x); plot(x, y 1, x, y 2, ’--’) %%%%%%%plot different lines plot(x, y 1, ’r’, x, y 2, ’k’) %%%%%change colour of line; %%%%%%%plot lines with different colours and styles plot(x, y 1, 'w') hold on plot(x, y 2, 'r', 'linestyle', '--') get(gcf) %%gives you all of the figure ‘handles’

Calculate the square root of x >>y 1=x^0. 5 or >>y 1=sqrt(x) Have to include a dot in order to multiply or divide element by element, x. ^0. 5 >>y 1=x. ^0. 5; y 2=x. ^1. 5; y 3=x. ^2; >>plot(x, y 1, x, y 2, ’—’, x, y 3, ’: ’)

To look at the importance of the ‘. ’ Go back to our matrix C 1=([2 1 4; 4 3 6; 6 8 2]) Now try to multiply C 1 by C 1 MESHGRID – spans a grid where points in a 2 D matrix are assigned X and Y positions (or X, Y, Z for 3 D) EG, we want to calculate the cosine Xby. Y on an area limited by X=+-3 and Y=+-3. Make the grid: >>[X, Y]=meshgrid(-3: . 1: 3, -3: . 1: 3) >>figure(1); Z=X. *Y; surf(X, Y, Z); colorbar %%%%%plots 3 D figure and adds a %%%%colorbar >>figure; Z=X. *Y; contourf(X, Y, Z); %%%%%2 D figure filled contours >>figure; Z=X. *Y; contour(X, Y, Z); %%%%%2 D figure not filled contours >> figure(2); Z=cos(X. *Y); surf(X, Y, Z);

Control flow of data…. . for i=1: 10, for i=10: 15 if i<5 counter=counter+1 a(i)=i*2; if i<12 a(counter)=i*2; elseif i==5 elseif i==12 a(i)=i; a(counter)=i; elseif i>12 a(i)=i/2; a(counter)=i*3; end end

>>while i>5; >> a(i)=i*2; >> i=i-1; >>end

Tasks 1) Create a vector X with values from 0 to 100 at intervals of 0. 1 a) What is the 20 th prime number in X b) How many prime numbers are less than 1 million c) Factorize 646969323 (i. e. find the sequence of prime numbers resulting in this number when multiplied together) 2) Create the vector X=[65 114 101 32 121 117 32 99 114 97 122 121 63 10] This is actually a message as each number represents a character or symbol. Matlab uses ASCII coding and every character is assigned a value between 0 and 255. a) Which sentence is represented by X? b) What is the ascii code for ‘enter’ c) Using findstr what is the position of r in the above text?

Code X=[0: 0. 1: 100]; t=primes(100); t(20) t 2=primes(1000000); size(t 2) f=factor(646969323); x=[65 114 101 32 121 117 32 99 114 97 122 121 63 10]; s=char(x); double('enter'); findstr(x, 'r')

Often need to work with 3 or 4 D matrices, such as global atmosphere data or model generated data. Here is an example where we create a 3 D space spanned by axes X, Y and Z. The first axis has values 0 to 360 with steps of 20. The second axis -90 to 90 with steps of 20 And finally an axis from 1990 to 1999. 99 with steps of 1/12: >>[X, Y, Z]=meshgrid(0: 20: 360, -90: 20: 90, 1990: 1/12: 1999. 99); Generate the hypothetical data……. Answer the following questions >>data=sin(X/1000). *cos(Y/100). *tanh(Z/2000);

a) What is the size of the matrix data? What is the minimum, maximum and mean value. In what position is the maximum? b) Assume that the first 2 dimensions represent longitude and latitude and the third is time and that the data are a scalar field eg global mean temps for a decade. How can you calculate the mean for each year? c) It is common to make anomalies of data (deviation from the climate mean) by removing the mean for every month. i. e for every january we remove the mean of all the January months in the data (hear we have 10). Jan. Anom 1999 = Janvalue 1999 - Janmean(1990 -1999); d) Plot the July anomaly for 1995 (help contour)

Answer Codes >> size(data) >>[val, ind]=max(data(: )) >>[a, b, c, ]=ind 2 sub(size(data), ind) ans = >>data 2=reshape(data, 10, 19, 12, 10); >>size(data 2) >>for y=1: 10; >> Yr_data=squeeze(data 2(: , : , y)); >> Yr_data 2(: , y)=mean(Yr_data, 3); >>end 10 19 120 >> min(data(: )) ans = 0 >> max(data(: )) ans = 0. 26694353105108 >> mean(data(: )) ans = 0. 11407671337753 >>for m=1: 12; >> mn_data=squeeze(data 2(: , m, : )); >> mn_climo(: , m)=mean(mn_data, 3); >>end >>for m=1: 12; >> for yr=1: 10 >> Mn_anom(: , m, yr)=squeeze(data 2(: , m, yr))squeeze(mn_climo(: , m)); >> end >>end

Reading ascii data w/o header Example ctd >> clear all >> load examplectd. dat Whos – examplectd in workspace which is 2138 x 5 Ignores spaces. But if data is separated by comma and NO SPACE, have to use function dlmread – data=dmlread(‘filename’, ’, ’); ALSO, if file has no extension have to include –ascii when loading…. >>load ‘filename’ –ascii Now, to split into variables: >>depth=examplectd(: , 1); temp=examplectd(: , 2); salinity=examplectd(: , 3); oxygen=examplectd(: , 4); fluro=examplectd(: , 5);

figure subplot(1, 4, 1); plot(temp, -depth); subplot(1, 4, 2); plot(salinity, -depth); subplot(1, 4, 3); plot(oxygen, -depth); subplot(1, 4, 4); plot(fluro, -depth); figure, plot(salinity, temp) title('salinity temp plot from BATS') xlabel('salinity') ylabel('temperature') legend

Ascii-data with boring header Use the function textread >> [cruise, decimal, lat, lon, press, depth, temp]=textread('examplectdtxt. d at', '%f %f', 'headerlines', 2); If csv >> [cruise, decimal, lat, lon, press, depth, temp]=textread('examplectdtxt. d at', '%f %*s %f %*s %f', 'headerlines', 2); The sign %* ensures that the data are not returned by the function. The last part of the command ‘headerlines’, 2 tells matlab to ignore the first 2 lines

Ascii-data with interesting header >>fid=fopen(‘examplectdtxt. dat’, ’r’); >>dummy=fread(fid, ’char’); >>fclose(fid); First we open file for reading (‘r’ is short for read). And a file identification number is given to fid. 2 nd line, we read all the data as ascii code (numbers between 0 and 127) And then on 3 rd line we close the file.

So we have now retrieved data and stored as numbers >> char(dummy(1: 1400)') Converts ascii code to characters. Note the ‘ transpose. With whos you can see that dummy is a 86069 x 1 matrix (column vector). So to print data to look like in a text editor, need to transpose to 1 x 86069.

Now we want to extract information from the header: year, lat and lon. Convert dummy into a row of characters >> c=char(dummy'); Then use findstr to find the positions of the variables we want >>temp_index=findstr('(C)', c); >>ctd=str 2 num(c(temp_index+3: end));

saving ascii-data The opposite function of load is save >>load examplectd. dat; >>save 'C: Documents and Settingsfialhoctd. asc' -ascii; if –ascii is not included, matlab will save as a. mat file, which is fine for use in matlab but not anything else. The new file is twice as large as the old one as matlab saves in a different format (look at two files). To separate data with column… >>dlmwrite('ctd 2. asc', examplectd, ', '); This is more compressed, but can be slow. A more readable form can be made using the fprintf function. Similar to textread

To separate data with column usinf fprintf… load 'examplectd. dat' >> fid=fopen('ctd 3. asc', 'wt') >> fprintf(fid, '%4 d %6. 3 f %6. 2 f %5. 3 fn', examplectd'); >> fclose(fid)

Plotting the CTD data from examplectd. dat >>load 'examplectd. dat' >> depth =examplectd(: , 1); >> temp =examplectd(: , 2); >> salinity=examplectd(: , 3); >> oxygen=examplectd(: , 4); >> fluoro =examplectd(: , 5); The figure window is at the top of the heirarchy! There are 62 parameters controlling the figure window. 2 commands used to change the settings: get and set >>figure, plot(temp, -depth) >>get(1) OR get(gcf) where gcf=get current figure. % lists figure parameters

The third parameter is Color = [0 0 0] telling us that the color of the figure is black >>get(gcf, ’color’); % change the color by setting parameter >>set(gcf, ’color’, [0. 2]) >>set(gcf, ’name’, ’ctd 1 temp’); %%%puts name on fig >>set(gcf, position’, [100 150 600 400]); %%can also use mouse!

axes that appear in the figure also have properties (~90), but now have to get all of the axis properties: >>get(gca) If you want red instead of black axis window: >>set(gca, ’color’, [1 0 0]); %%%%%%%other parameters to adjust: >>set(gca, ’fontangle’, ’italic’) >>set(gca, ’fontsize’, 15) >>set(gca, ’linewidth’, 2) >>set(gca, ’position’, [. 2. 2. 6. 6]) >>set(gca, ’xcolor’, [0 0 1]) >>set(gca, ’xlim’, [1 24]) >>set(gca, ’xtick’, [1: 5: 24]) >>set(gca, 'xticklabel', ['01 C'; '06 C'; '11 C'; '16 C'; '21 C'])

can do many at the same time: >> set(gca, 'fontangle', 'italic', 'fontsize', 15, 'linewidth', 2, 'xgrid', 'on') Adding text: >> title('ctd temperature plot') >> xlabel('temperature') >> ylabel('Depth (m)') >> gtext('18 C water') %%%%puts the text where you click on fig also the functions: grid on and grid off hold on and hold off

Plotting the Current. Meter data from adcp. txt >> load 'ADCP. txt' -ascii >> direction=ADCP(: , 1); >> Vel. E=ADCP(: , 2); >> Vel. N=ADCP(: , 3); >> Vel. U=ADCP(: , 4); >> figure, plot(Vel. U) >> figure, plot(Vel. E) >> figure, plot(Vel. N) other plotting options available: polar, rose and feather in matlab the direction is needed in radians, and often the direction is given relative to north…… to find the direction in radians: >>theta=(90 -direction)*pi/180;

Plotting the SST data from SST_IO. mat >> load('SST_IO') >> whos >> climo. SST=squeeze(mean(SST)); >> for mn=1: 12; SSTanom(mn, : )=squeeze(SST(mn, : ))-climo. SST; end %%%now we want to plot some contoured figures >>temprange=[20: 0. 5: 30] %%%defines the contours we’ll plot >>for mn=1: 12 >> subplot(3, 4, mn) >> contourf(lon, lat, squeeze(SST(mn, : )), temprange) >> caxis([20 30]); %%%%%%defines the limits of the colorbar >>end >>colorbar

Or maybe we just want to plot a figure showing the annual mean and we want to interpolate >> figure >> pcolor(lon, lat, climo. SST); shading interp >> colorbar %%%%%then we want to plot contourf lines on top for 25 C >> hold on >> [cs, h]=contour(lon, lat, climo. SST, [27 28 29], 'w'); %%%[cs, h] is a way %%of giving the contour lines a ‘handle’ in order to label them >>clabel(cs, ’manual’);

%%now lets look at the anomalies >>temprange=[-10: 0. 2: 10] %%%defines the contours we’ll plot >>for mn=1: 12 subplot(3, 4, mn) contourf(lon, lat, squeeze(SSTanom(mn, : )), temprange) caxis([-2 2]); %%%%%%defines the limits of the colorbar end >>colorbar look at times series over IO >> anom 1=squeeze(nanmean(SSTanom, 2)); %%%use nanmean as %%there are Na. N values as a landmask >> anom. IO=squeeze(nanmean(anom 1, 2)); >> figure, plot(anom. IO)

%%%%%%Do some simple correlations, so create a similar time-series this is just a nonsensical timeseries we will create and call atlantic ocean. Obviously this is just a timeseries we create from the indian ocean data but it saves downloading further data! (a bit nonsensical but then you get the correlation programs!!!) >> AO_anom=tan(anom. IO); >> figure, plot(AO_anom), hold on, plot(anom. IO, ’: ’)

y 1=AO_anom; y 2=anom. IO; %%%%input variables clear brol 1 brol 2 r 1 r 2 tau eodf r_95 prob t n r; %%%%%clear specific variables (if run this before) r=corrcoef(y 1, y 2); %%%%%%%calcluate correlation [c, lags]=xcorr(y 1, y 2, 'coeff'); %%%%%calculate correlation at different time lags between data if (abs(min(c))>=max(c)); %%%%find if the max correlation is positive or negative r=min(c), [x, y]=min(c); else; r=max(c), [x, y]=max(c); end n=length(y 1); %%%%%%%% k=n/2; %%%%%%%% for i=1: k; brol 1=corrcoef(y 1(1: n-i), y 1(i+1: n)); brol 2=corrcoef(y 2(1: n-i), y 2(i+1: n)); r 1(i)=brol 1(1, 2); r 2(i)=brol 2(1, 2); end; tau=[1+2*sum(r 1. *r 2)]; eodf=round(n/tau); r_95=1. 96/sqrt(eodf-2) %%%%%%correlation at. 95 confidence r_99=2. 33/sqrt(eodf-2) %%%%%%correlation at. 99 confidence t= abs(r/sqrt((1 -r*r)/(eodf-2))); df=eodf-2; %prob_directional=tcdf(t, df) %prob_non_directional=1 -(1 -prob_directional)*2 r %correlation, is it significant at the. 95 or. 99 level? lags(y) %%%% time lag (based on the timesteps of data set, if -ve then y 2 leads y 1

%%%%%%%help m_map – built in toolbox which allows you to plot on different projections (useful if looking at large areas). Also has built in topography and coastlines etc