MATLAB Cheatsheet

% row vector
v = [1 2 3 4 5];

% column vector (semicolon = new row)
c = [1; 2; 3];

% 3x3 matrix
A = [1 2 3; 4 5 6; 7 8 9];

% special matrices
Z = zeros(3);          % 3x3 zeros
O = ones(2, 4);        % 2x4 ones
I = eye(3);            % 3x3 identity
R = rand(2, 3);        % 2x3 uniform random
N = randn(4);          % 4x4 normal random
M = magic(4);          % 4x4 magic square
D = diag([1 2 3]);     % diagonal matrix

% linear spacing
l = linspace(0, 1, 5); % [0 0.25 0.5 0.75 1]
r = 0:0.5:2;           % [0 0.5 1 1.5 2]

A = magic(4);          % 4x4 magic square

% single element (1-indexed!)
A(2, 3)                % row 2, col 3

% entire row or column
A(1, :)                % first row (all columns)
A(:, 2)                % second column (all rows)

% submatrix
A(1:2, 2:4)            % rows 1-2, cols 2-4

% linear indexing (column-major)
A(5)                   % 5th element counting down columns

% end keyword
A(end, :)              % last row
A(:, end-1)            % second-to-last column

% logical indexing
A(A > 10)              % all elements > 10 (as vector)
A(A > 10) = 0;         % set large elements to 0

% find indices
[r, c] = find(A == 16) % row and col of value 16

A = [1 2; 3 4];
B = [5 6; 7 8];

% matrix operations
C = A + B;             % addition
D = A - B;             % subtraction
E = A * B;             % matrix multiplication
F = A';                % transpose
G = A^2;               % matrix power (A * A)

% element-wise operations (use the dot!)
H = A .* B;            % element-wise multiply
I = A ./ B;            % element-wise divide
J = A .^ 2;            % element-wise power
K = 2 * A;             % scalar multiply

% matrix functions
det(A)                 % determinant
inv(A)                 % inverse
pinv(A)                % pseudo-inverse
rank(A)                % rank
trace(A)               % trace (sum of diagonal)

% solving linear systems Ax = b
b = [5; 11];
x = A \ b;              % left division: solves A*x = b
x = inv(A) * b;        % equivalent but slower

A = [1 2; 3 4];
B = [5 6; 7 8];

% concatenation
C = [A B];             % horizontal: [1 2 5 6; 3 4 7 8]
D = [A; B];            % vertical: 4x4

% reshape
E = reshape(1:12, 3, 4);  % 3x4 matrix from 1..12

% flip and rotate
F = flipud(A);         % flip up-down
G = fliplr(A);         % flip left-right
H = rot90(A);          % rotate 90 degrees

% repmat (tile)
I = repmat([1 2], 2, 3);  % repeat [1 2] in a 2x3 grid

% size and length
[m, n] = size(A);      % m=2, n=2
len = length(A);       % max dimension = 2
num = numel(A);        % total elements = 4

% remove rows/cols
A(2, :) = [];          % delete row 2
A(:, 1) = [];          % delete column 1

% math functions operate element-wise
x = 0:0.1:2*pi;
y = sin(x);            % sine of each element
z = exp(x);            % exponential
w = sqrt(x);           % square root

% rounding
round(3.7)             % 4
floor(3.7)             % 3
ceil(3.2)              % 4
fix(-3.7)              % -3 (toward zero)

% sums and products along dimensions
A = magic(3);
sum(A)                 % sum of each column (row vector)
sum(A, 2)              % sum of each row (column vector)
sum(A(:))              % sum of ALL elements
prod(A)                % product of each column
cumsum(A)              % cumulative sum

% min, max, sort
[v, i] = max(A(:))     % max value and its index
sort(A, 'descend')     % sort each column descending

% AVOID loops — vectorize!
% BAD: for i=1:1000, y(i)=sin(i); end
% GOOD:
i = 1:1000;
y = sin(i);            % fast, vectorized

score = 85;

if score >= 90
    grade = 'A';
elseif score >= 80
    grade = 'B';
elseif score >= 70
    grade = 'C';
else
    grade = 'F';
end

disp(grade);  % B

% logical operators: && (and), || (or), ~ (not)
if x > 0 && x < 10
    disp('in range');
end

% compare strings with strcmp
if strcmp(name, 'Alice')
    disp('hi Alice');
end

% for loop over a range
for i = 1:5
    fprintf('i = %d\n', i);
end

% iterate over a vector
v = [10 20 30];
for val = v
    disp(val);
end

% nested loop over a matrix
A = zeros(3, 3);
for r = 1:3
    for c = 1:3
        A(r, c) = r * c;
    end
end

% while loop
n = 10;
while n > 1
    n = n / 2;
    fprintf('%.2f\n', n);
end

% preallocate for speed (IMPORTANT!)
result = zeros(1, 1000);
for i = 1:1000
    result(i) = i^2;
end

% switch statement
method = 'linear';
switch method
    case 'linear'
        disp('using linear');
    case 'cubic'
        disp('using cubic');
    case {'nearest', 'spline'}
        disp('using nearest or spline');
    otherwise
        disp('unknown method');
end

% break and continue
for i = 1:10
    if i == 5
        break        % exit the loop entirely
    end
    if mod(i, 2) == 0
        continue     % skip to next iteration
    end
    disp(i);         % prints 1 3
end

% try-catch for error handling
try
    x = 1 / 0;       % may error
catch ME
    fprintf('Error: %s\n', ME.message);
end

A = [1 2 3 4 5 6 7 8 9 10];

% comparison operators
A > 5                  % logical array
A == 5
A ~= 5                 % not equal (NOT !=)
A >= 3 & A <= 7        % element-wise AND
A < 3 | A > 7          % element-wise OR

% logical indexing (powerful!)
A(A > 5)               % [6 7 8 9 10]
A(A > 5) = 0;          % set elements > 5 to 0
A(mod(A, 2) == 0)      % even numbers

% find indices
idx = find(A > 5);     % indices of elements > 5
[r, c] = find(A > 5);  % row and col (for matrices)

% any and all
any(A > 5)             % true if ANY element > 5
all(A > 0)             % true if ALL elements > 0
any(A > 5, 2)          % per row

% is functions
isnan(x); isinf(x); isfinite(x);
isnumeric(x); ischar(x); iscell(x);

% string creation
s1 = 'Hello';          % char array
s2 = "World";          % string scalar (R2017+)

% concatenation
greeting = [s1 ', ' s2 '!'];  % char array concat
full = s1 + " " + s2;          % string concat

% formatting
name = 'Alice';
age = 30;
fprintf('Name: %s, Age: %d\n', name, age);
str = sprintf('Pi is %.2f', pi);  % returns string

% common string functions
upper('hello')         % HELLO
lower('WORLD')         % world
strlength("hello")     % 5 (string)
length('hello')        % 5 (char array)
strcmp('a', 'a')       % true
strfind('hello', 'll') % 3 (index)
strrep('cat', 'c', 'b') % bat
strsplit('a,b,c', ',') % {'a','b','c'}
strjoin({'a','b'}, '-') % a-b
num2str(42)             % '42'

% function in a file named 'add.m'
function result = add(a, b)
    % ADD returns the sum of a and b
    result = a + b;
end

% with multiple inputs and outputs
function [mn, mx, avg] = stats(v)
    mn = min(v);
    mx = max(v);
    avg = mean(v);
end

% calling functions
s = add(3, 4);              % 7
[lo, hi, mu] = stats([1 2 3 4 5]);

% capture only first output
minimum = stats([1 2 3]);   % gets mn only

% nargin and nargout (number of args)
function y = power(x, n)
    if nargin < 2
        n = 2;              % default value
    end
    y = x .^ n;
end

% anonymous function (one-liner, no file needed)
square = @(x) x .^ 2;
square(5)              % 25
square([1 2 3])        % [1 4 9]

% with multiple inputs
add = @(a, b) a + b;
add(3, 4)              % 7

% capture variables from workspace
c = 10;
addc = @(x) x + c;
addc(5)                % 15

% function handle to built-in
f = @sin;
f(pi/2)                % 1

% pass functions as arguments
result = integral(@(x) x.^2, 0, 1);  % integrate x^2 from 0 to 1
result = fzero(@(x) x^2 - 2, 1);     % find root near 1

% array of function handles
funs = {@sin, @cos, @tan};
funs{1}(0)             % 0

% A script is just a sequence of commands in a .m file
% It shares the base workspace

% script: analyze_data.m
data = load('data.mat');
cleaned = data.values(data.values > 0);
mean_val = mean(cleaned);
fprintf('Mean: %.2f\n', mean_val);
plot(cleaned);
title('Cleaned Data');

% sections (cells) with %%
%% Initialize
x = linspace(0, 2*pi, 100);

%% Plot
plot(x, sin(x));

%% Analyze
disp(mean(sin(x)));

% Run a section: Ctrl+Enter (in editor)
% Live Scripts (.mlx): rich text, inline plots, equations

% local functions in a script (must be at the end)
% main_script.m
x = 1:10;
y = process(x);
disp(y);

function r = process(v)
    r = normalize(scale(v));  % calls another local function
end

function s = scale(v)
    s = v * 10;
end

function n = normalize(v)
    n = v / max(v);
end

% nested function (inside another function)
function outer(x)
    y = 0;
    function inner()
        y = y + x;    % can access outer's variables
    end
    inner();
    disp(y);
end

% base workspace variables
x = 10;

% functions have their OWN workspace
function y = myfunc()
    % x is NOT accessible here
    y = 5;
end

% pass by value (modifications don't affect caller)
function y = modify(v)
    v = v * 2;        % local copy
    y = v;
end
a = 5;
b = modify(a);        % b=10, a still 5

% global variables (use sparingly!)
global COUNTER
COUNTER = 0;
function increment()
    global COUNTER
    COUNTER = COUNTER + 1;
end

% persistent variables (like static in C)
function counter()
    persistent count
    if isempty(count)
        count = 0;
    end
    count = count + 1;
    disp(count);
end

% cell arrays hold mixed data types
c = {'Alice', 30, [1 2 3], magic(2)};

% access with {} (content) vs () (cell)
c{1}                  % 'Alice' (the string itself)
c(1)                  % {'Alice'} (a 1x1 cell)
c{3}(2)               % 2 (second element of third cell)

% modify and append
c{2} = 31;
c{end+1} = 'new';     % append
c{5} = 'skip';        % auto-fills gaps with []

% multi-dimensional cell array
c2 = {'a', 'b'; 'c', 'd'};  % 2x2 cell

% iterate over cells
for i = 1:numel(c)
    disp(c{i});
end

% convert: cell2mat, mat2cell, num2cell
nums = num2cell(1:5);  % {1, 2, 3, 4, 5}

% struct: named fields with any data
user.name = 'Alice';
user.age = 30;
user.scores = [90 85 88];
disp(user.name);      % Alice

% struct constructor
s = struct('name', 'Bob', 'age', 25, 'active', true);

% array of structs
users(1) = struct('name', 'Alice', 'age', 30);
users(2) = struct('name', 'Bob', 'age', 25);
disp(users(2).name);  % Bob

% table (like a dataframe) — great for tabular data
names = {'Alice'; 'Bob'; 'Carol'};
ages = [30; 25; 28];
scores = [90; 85; 88];
T = table(names, ages, scores);

% access table columns
T.ages               % column as array
T(:, 2)              % column as table
T(1:2, :)            % first 2 rows
T.age > 26           % logical indexing
T(T.age > 26, :)     % rows where age > 26

% summary statistics
summary(T)

% string array (R2017+) — recommended for text
names = ["Alice", "Bob", "Carol"];
names(1)              % "Alice"
names'                % column string array
strlength(names)      % [5 3 5]

% char array (legacy)
c = 'Hello World';
c(1:5)                % 'Hello'
size(c)               % [1 11]

% cell array of char vectors (legacy mixed-length)
old = {'Alice', 'Bob', 'Carol'};

% convert between types
s = string('hello');          % char -> string
c = char("hello");            % string -> char
cellarr = cellstr(names);     % string array -> cell of char

% combine and split
full = "Alice" + " " + "Smith";  % "Alice Smith"
parts = split("a,b,c", ",");     % ["a", "b", "c"]
joined = join(["a" "b" "c"], "-"); % "a-b-c"

% pattern matching
matches(names, 'A*')  % logical: [true false false]
contains(names, 'li') % [true false false]
replace("hello", 'l', 'L')  % "heLLo"

% map (dictionary / hash table)
m = containers.Map;
m('name') = 'Alice';
m('age') = 30;
m('scores') = [90 85];

% access
m('name')             % 'Alice'
isKey(m, 'name')      % true
keys(m)               % {'age', 'name', 'scores'}
values(m)             % {30, 'Alice', [90 85]}
remove(m, 'age');

% create with initial values
m2 = containers.Map({'a','b','c'}, {1, 2, 3});

% set operations
A = [1 2 3 4 5];
B = [4 5 6 7 8];
union(A, B)           % [1 2 3 4 5 6 7 8]
intersect(A, B)       % [4 5]
setdiff(A, B)         % [1 2 3]
setxor(A, B)          % [1 2 3 6 7 8]
ismember(3, A)        % true
ismember([1 6], A)    % [true false]

% categorical data (efficient for repeated strings)
colors = categorical({'red', 'blue', 'red', 'green', 'blue'});
categories(colors)    % list of categories
summary(colors)       % count per category
ord = categorical({'low','med','high'}, ...
     {'low','med','high'}, 'Ordinal', true);
ord(1) < ord(2)       % true (ordered)

% datetime
t = datetime('now')   % current date and time
t = datetime(2024, 1, 15)
t.Format = 'yyyy-MM-dd HH:mm:ss';

% duration
d = duration(2, 30, 0);  % 2 hours 30 min
d = hours(2) + minutes(30);

% date arithmetic
t2 = t + days(7);     % add a week
t2 - t                % duration of 7 days

% date strings and numbers
s = datestr(t, 'yyyy-mm-dd');  % to string (legacy)
n = datenum(t);                 % to serial number
t = datetime(n, 'ConvertFrom', 'datenum');

% generate date range
dates = datetime(2024,1,1) + days(0:6);  % one week

x = linspace(0, 2*pi, 100);

% line plot
plot(x, sin(x), 'b-', 'LineWidth', 2);
hold on;
plot(x, cos(x), 'r--', 'LineWidth', 2);
hold off;

% labels and legend
xlabel('x (radians)');
ylabel('amplitude');
title('Sine and Cosine');
legend('sin(x)', 'cos(x)', 'Location', 'best');
grid on;

% line styles: b- blue solid, r-- red dashed
% g: green dotted, m-. magenta dash-dot
% markers: o circle, * star, s square, + plus

% axis control
axis([0 2*pi -1.5 1.5]);  % [xmin xmax ymin ymax]
xlim([0 2*pi]);
ylim([-1.5 1.5]);

% save figure
saveas(gcf, 'plot.png');
print('-dpng', '-r300', 'plot.png');

% subplot(rows, cols, index)
figure;
subplot(2, 2, 1);
plot(x, sin(x));
title('Sine');

subplot(2, 2, 2);
plot(x, cos(x));
title('Cosine');

subplot(2, 2, 3);
plot(x, tan(x));
title('Tangent');
ylim([-5 5]);

subplot(2, 2, 4);
plot(x, exp(-x));
title('Decay');

% tiledlayout (newer, cleaner)
figure;
t = tiledlayout(2, 2);
nexttile; plot(x, sin(x)); title('Sine');
nexttile; plot(x, cos(x)); title('Cosine');
nexttile; plot(x, tan(x)); title('Tangent');
nexttile; plot(x, exp(-x)); title('Decay');
title(t, 'Trig Functions');  % overall title

% multiple y-axes
yyaxis left;  plot(x, sin(x));
yyaxis right; plot(x, 100*cos(x));

% bar chart
bar([1 2 3], [10 15 7]);
bar([1 2 3; 4 5 6], 'grouped');  % grouped bars
bar([1 2 3; 4 5 6], 'stacked');  % stacked bars

% histogram
data = randn(1000, 1);
histogram(data, 30);    % 30 bins
histogram(data, 'Normalization', 'pdf');

% scatter plot
x = rand(100, 1);
y = 2*x + randn(100, 1)*0.1;
scatter(x, y, 50, 'filled');
% colored scatter
scatter(x, y, 50, y, 'filled');  % color by y value
colorbar;

% pie chart
pie([30 20 50], {'A', 'B', 'C'});

% area plot
area(1:5, [1 3 2 4 5]);

% stem (discrete)
stem(0:5, [1 4 9 16 25 36]);

% boxplot (requires Statistics toolbox)
boxplot(randn(100, 3));

% 3D line plot
t = 0:0.1:10*pi;
plot3(sin(t), cos(t), t);
xlabel('x'); ylabel('y'); zlabel('z');

% meshgrid for surface plots
[X, Y] = meshgrid(-2:0.1:2);
Z = X .* exp(-X.^2 - Y.^2);

% surface plots
figure;
subplot(1, 3, 1);
surf(X, Y, Z);     % surface
title('surf');

subplot(1, 3, 2);
mesh(X, Y, Z);     % wireframe
title('mesh');

subplot(1, 3, 3);
contour(X, Y, Z, 20);  % contour lines
title('contour');

% colored surface with colorbar
figure;
surf(X, Y, Z, 'EdgeColor', 'none');
colorbar;
colormap jet;
shading interp;     % smooth colors

% view angle
view(45, 30);       % azimuth, elevation
axis equal;

% customize line appearance
plot(x, y, ...
    'Color', [0.2 0.6 0.8], ...     % RGB 0-1
    'LineWidth', 2, ...
    'Marker', 'o', ...
    'MarkerSize', 8, ...
    'MarkerFaceColor', 'r');

% annotations
text(1, 0.5, 'important point');
annotation('arrow', [0.3 0.5], [0.3 0.5]);

% figure and axes properties
set(gcf, 'Position', [100 100 800 600]);
set(gca, 'FontSize', 14, 'FontName', 'Arial');
set(gca, 'XScale', 'log');  % log scale

% latex in labels
title('Function: \\alpha + \\beta^2', 'Interpreter', 'latex');

% export high quality
exportgraphics(gcf, 'plot.pdf', 'ContentType', 'vector');
exportgraphics(gcf, 'plot.png', 'Resolution', 300);

% legend customization
legend('show', 'Location', 'northwest', ...
       'FontSize', 12, 'Box', 'off');

% multiple figures
f1 = figure; plot(x, sin(x));
f2 = figure; plot(x, cos(x));
figure(f1);  % switch back to f1

data = [1 2 3 4 5 6 7 8 9 10];

% central tendency
mean(data)            % 5.5
median(data)          % 5.5
mode(data)            % 1 (first mode)
geomean(data)         % geometric mean (Stats toolbox)
harmmean(data)        % harmonic mean

% dispersion
std(data)             % standard deviation
var(data)             % variance
range(data)           % max - min
iqr(data)             % interquartile range
mad(data)             % mean absolute deviation

% extremes and quantiles
min(data); max(data);
quantile(data, 0.25)  % 25th percentile
prctile(data, 95)     % 95th percentile

% along dimensions
A = randn(100, 3);
mean(A)               % mean of each column
mean(A, 2)            % mean of each row
mean(A, 'all')        % overall mean

% ignoring NaNs
data(5) = NaN;
mean(data, 'omitnan') % skip NaN values
nanmean(data)         % legacy equivalent

% polynomial fit
x = 0:0.5:10;
y = 2*x + 1 + randn(size(x))*2;  % noisy linear data

% linear least squares (degree 1)
p = polyfit(x, y, 1);   % [slope, intercept]
yfit = polyval(p, x);
plot(x, y, 'o', x, yfit, '-');

% higher-order polynomial
p3 = polyfit(x, y, 3);
yfit3 = polyval(p3, x);

% fit with fittype (Curve Fitting Toolbox)
f = fit(x', y', 'poly2');  % quadratic
plot(f, x, y);

% custom model fit
% fit y = a*exp(b*x) using lsqcurvefit
model = @(p, x) p(1)*exp(p(2)*x);
params = lsqcurvefit(model, [1 0.1], x, y);

% correlation
r = corrcoef(x, y);     % correlation matrix
r = corr(x, y);         % correlation coefficient

% linear regression with stats
X = [ones(size(x)) x'];  % design matrix
[b, bint, r, rint, stats] = regress(y', X);
% stats = [R^2, F, p, error_var]

% 1D interpolation
x = 0:5;
y = [0 1 4 9 16 25];
xi = 0:0.5:5;

yi = interp1(x, y, xi, 'linear');   % linear (default)
yi = interp1(x, y, xi, 'spline');   % cubic spline
yi = interp1(x, y, xi, 'pchip');    % shape-preserving
yi = interp1(x, y, xi, 'nearest');  % nearest neighbor

% extrapolation
yi = interp1(x, y, xi, 'linear', 'extrap');

% 2D interpolation
[X, Y] = meshgrid(0:2:10);
Z = X.^2 + Y.^2;
[Xi, Yi] = meshgrid(0:0.5:10);
Zi = interp2(X, Y, Z, Xi, Yi, 'spline');

% resample signal to new rate
y_resampled = resample(y, 3, 2);  % 3/2 times the rate

% griddedInterpolant (efficient for repeated queries)
F = griddedInterpolant(x, y, 'spline');
y1 = F(2.5);
y2 = F([1.5 3.5 4.5]);

% generate a signal with two frequencies
fs = 1000;                    % sampling rate (Hz)
t = 0:1/fs:1-1/fs;            % 1 second of time
x = 2*sin(2*pi*50*t) + 1*sin(2*pi*120*t);

% compute FFT
N = length(x);
Y = fft(x);
P2 = abs(Y/N);                % two-sided spectrum
P1 = P2(1:N/2+1);             % one-sided spectrum
P1(2:end-1) = 2*P1(2:end-1);
f = fs*(0:(N/2))/N;

% plot
plot(f, P1);
xlabel('Frequency (Hz)');
ylabel('Amplitude');
title('Single-Sided Spectrum');

% filtering (simple low-pass)
cutoff = 100;  % Hz
Y_filt = Y;
Y_filt(f > cutoff) = 0;       % zero out high freqs (simplified)
x_filtered = real(ifft(Y_filt));

% spectrogram
spectrogram(x, 128, 120, 128, fs, 'yaxis');

% find root of f(x) = 0
f = @(x) x^2 - 2;
x_root = fzero(f, 1);         % ~1.4142 (sqrt(2))
x_root = fsolve(f, 1);        % alternative

% minimize a function
fun = @(x) (x-3)^2 + 1;
[x_min, f_min] = fminbnd(fun, 0, 5);  % bounded 1D
% x_min = 3, f_min = 1

% multivariate minimization
fun2 = @(x) x(1)^2 + x(2)^2;  % sphere function
x0 = [1, 1];
[x_opt, f_opt] = fminsearch(fun2, x0);
% x_opt ≈ [0, 0]

% constrained optimization (Optimization Toolbox)
% minimize fun subject to A*x <= b
A = []; b = []; Aeq = []; beq = [];
lb = [0 0]; ub = [];
[x, fval] = fmincon(fun2, x0, A, b, Aeq, beq, lb, ub);

% linear programming
% minimize f'*x subject to A*x <= b
f = [-1; -1];                 % maximize x+y
A = [1 1; -1 0; 0 -1];
b = [2; 0; 0];
[x, fval] = linprog(f, A, b);

% save variables to a .mat file
x = 1:10;
y = sin(x);
save('data.mat');             % save ALL variables
save('data.mat', 'x', 'y');   % save specific variables
save('data.mat', '-append');  % add more variables

% load variables
load('data.mat');             % load all into workspace
loaded = load('data.mat');    % load into a struct
loaded.x                      % access via struct

% save with compression
save('data.mat', 'x', 'y', '-v7.3');  % large files (>2GB)
save('data.mat', '-v7');              % compressed (default)

% save specific format
save('data.txt', 'x', '-ascii');      % plain text
save('data.csv', 'x', '-ascii', '-double');

% clear variables
clear x y;                    % remove specific
clear all;                    % remove everything
clearvars -except x;          % keep only x

% readtable (recommended for tabular data)
T = readtable('data.csv');
T = readtable('data.xlsx', 'Sheet', 'Sheet1');
T = readtable('data.csv', 'Delimiter', ',');

% access columns by name
T.Var1; T.age;

% readmatrix (numeric data)
M = readmatrix('numbers.csv');
M = readmatrix('data.txt', 'Delimiter', '\t');

% readcell (mixed data)
C = readcell('mixed.csv');

% low-level file I/O
fid = fopen('data.txt', 'r');
while ~feof(fid)
    line = fgetl(fid);
    disp(line);
end
fclose(fid);

% read formatted data
fid = fopen('data.txt', 'r');
data = fscanf(fid, '%f %f', [2, inf]);
fclose(fid);

% write text
fid = fopen('output.txt', 'w');
fprintf(fid, 'Result: %.2f\n', 3.14159);
fclose(fid);

% read Excel files
T = readtable('data.xlsx');
T = readtable('data.xlsx', 'Sheet', 'Sales', 'Range', 'A1:D100');

% numeric data only
M = readmatrix('data.xlsx', 'Sheet', 1);

% write to Excel
writetable(T, 'output.xlsx');
writetable(T, 'output.xlsx', 'Sheet', 'Results');

% write matrix
writematrix(M, 'output.xlsx', 'Sheet', 1);

% cell array (mixed types)
C = {'Name','Age'; 'Alice',30; 'Bob',25};
writecell(C, 'people.xlsx');

% append to existing sheet
writetable(T, 'data.xlsx', 'Sheet', 1, 'WriteMode', 'append');

% read specific range
T = readtable('data.xlsx', 'Range', 'B2:D10');

% detect import options (customize parsing)
opts = detectImportOptions('data.csv');
opts = setvartype(opts, 'age', 'int32');
T = readtable('data.csv', opts);

% current directory
pwd
cd('C:\\Users\\data')

% list files
files = dir();              % struct array of files
files = dir('*.mat');       % filter by pattern
files = dir('subdir/');     % list a subdirectory

% access file info
for i = 1:length(files)
    if ~files(i).isdir
        fprintf('%s - %d bytes\n', files(i).name, files(i).bytes);
    end
end

% path operations
fullfile('data', '2024', 'file.csv')  % data/2024/file.csv
[filepath, name, ext] = fileparts('C:\\data\\test.csv')
% filepath = 'C:\data', name = 'test', ext = '.csv'

% create and remove directories
mkdir('output');
rmdir('output', 's');      % remove recursively

% file existence
if exist('data.mat', 'file')
    load('data.mat');
end

% add folder to path
addpath('C:\\myfunctions');
savepath();                 % save for future sessions

% read and display an image
img = imread('photo.jpg');
imshow(img);
[m, n, c] = size(img);     % dimensions, c=3 for RGB

% write image
imwrite(img, 'output.png');
imwrite(img, 'output.jpg', 'Quality', 90);

% convert to grayscale
if size(img, 3) == 3
    gray = rgb2gray(img);
    imshow(gray);
end

% image processing
img2 = imresize(img, 0.5);          % scale 50%
img3 = imrotate(img, 45);           % rotate 45 degrees
img4 = imadjust(img, [0.2 0.8], []); % contrast adjust

% audio
[y, fs] = audioread('song.wav');    % y = samples, fs = sample rate
sound(y, fs);                        % play audio
audiowrite('output.wav', y, fs);

% record audio
recObj = audiorecorder(44100, 16, 1);
recordblocking(recObj, 3);           % record 3 seconds
y = getaudiodata(recObj);

% create symbolic variables (Symbolic Math Toolbox)
syms x y z

% symbolic expressions
f = x^2 + 2*x + 1;
g = sin(x)^2 + cos(x)^2;

% simplify
simplify(f)             % (x + 1)^2
simplify(g)             % 1
expand((x+1)^3)         % x^3 + 3*x^2 + 3*x + 1
factor(x^2 - 1)         % (x - 1)*(x + 1)
collect(x^2 + 2*x + x^2)% 2*x^2 + 2*x

% substitute values
subs(f, x, 3)           % 16 (f at x=3)
subs(f, x, y)           % y^2 + 2*y + 1

% pretty print
pretty(f)
disp(f)

% convert between numeric and symbolic
double(sym('1/3'))      % 0.3333
sym(0.5)                % 1/2
vpa(pi, 50)             % 50-digit pi: 3.1415926535897...

syms x

% differentiation
f = x^3 + 2*x^2 + x;
df = diff(f)            % 3*x^2 + 4*x + 1
d2f = diff(f, 2)        % second derivative: 6*x + 4
diff(f, x, 3)           % third derivative: 6

% partial derivatives
syms x y
g = x^2 * y + sin(x*y);
diff(g, x)              % 2*x*y + y*cos(x*y)
diff(g, y)              % x^2 + x*cos(x*y)

% integration
int(x^2)                % x^3/3 (indefinite)
int(x^2, x, 0, 1)       % 1/3 (definite, 0 to 1)
int(sin(x), x, 0, pi)   % 2

% limits
limit(sin(x)/x, x, 0)   % 1
limit(1/x, x, inf)      % 0
limit(1/x, x, 0, 'left')  % -inf

% Taylor series
taylor(exp(x), x, 0, 'Order', 6)
% 1 + x + x^2/2 + x^3/6 + x^4/24 + x^5/120

syms x y

% solve algebraic equations
solve(x^2 - 4 == 0, x)        % [2; -2]
solve(x^2 + 1 == 0, x)        % [i; -i] (complex roots)

% solve a system
sol = solve([x + y == 5, x - y == 1], [x, y]);
sol.x                        % 3
sol.y                        % 2

% solve with parameters
syms a b c
r = solve(a*x^2 + b*x + c == 0, x);
% r = -(b + (b^2 - 4*a*c)^(1/2))/(2*a), ...

% numerical solve (no closed form)
vpasolve(x^5 - 3*x + 1 == 0, x)

% differential equations
syms y(t) t
ode = diff(y, t) == -k*y;     % dy/dt = -k*y
cond = y(0) == y0;
ySol(t) = dsolve(ode, cond);  % y0*exp(-k*t)

% system of ODEs
syms x(t) y(t)
eq1 = diff(x,t) == y;
eq2 = diff(y,t) == -x;
sol = dsolve([eq1, eq2]);

syms t s w

% Laplace transform
f = exp(-a*t);
F = laplace(f)          % 1/(a + s)
% default: f(t) -> F(s)

% inverse Laplace
f2 = ilaplace(1/(s+1))  % exp(-t)

% Fourier transform
g = exp(-t^2);
G = fourier(g)          % pi^(1/2)*exp(-w^2/4)

% inverse Fourier
g2 = ifourier(G)

% Z-transform (for discrete signals)
syms n z
h = a^n;
H = ztrans(h)           % -z/(a - z) (for |z| > |a|)
h2 = iztrans(z/(z-a))   % a^n

% apply transform to solve ODE
% y'' + y = 0, y(0)=0, y'(0)=1
syms y(t) Y(s)
ode = diff(y, t, 2) + y == 0;
% take Laplace of both sides, solve for Y(s), then invert

% solve dy/dt = f(t, y) numerically with ode45

% define the ODE as a function
f = @(t, y) -2*y;       % dy/dt = -2y (exponential decay)

% solve on [0, 5] with y(0) = 1
[t, y] = ode45(f, [0 5], 1);
plot(t, y);
xlabel('t'); ylabel('y');

% system of ODEs (Lotka-Volterra predator-prey)
% dx/dt = 1.5x - xy,  dy/dt = -0.7y + 0.1xy
ode_sys = @(t, v) [1.5*v(1) - v(1)*v(2);
                   -0.7*v(2) + 0.1*v(1)*v(2)];
[t, V] = ode45(ode_sys, [0 20], [10; 5]);
plot(t, V(:,1), 'b', t, V(:,2), 'r');
legend('Prey', 'Predator');

% stiff ODE solver (for fast dynamics)
[t, y] = ode15s(f, [0 5], 1);

% set options
opts = odeset('RelTol', 1e-6, 'AbsTol', 1e-9);
[t, y] = ode45(f, [0 5], 1, opts);

% pass parameters via anonymous function
a = 0.5;
f = @(t, y, a) -a*y;
[t, y] = ode45(@(t,y) f(t,y,a), [0 10], 1);

% solve Ax = b
A = [3 2 -1; 2 -2 4; -1 0.5 -1];
b = [1; -2; 0];
x = A \ b;             % backslash: preferred solver
x = inv(A) * b;         % explicit inverse (slower, less stable)
x = linsolve(A, b);     % optimized linear solver

% LU, QR, Cholesky decompositions
[L, U, P] = lu(A);      % PA = LU
[Q, R] = qr(A);         % QR factorization
R = chol(A'*A);         % Cholesky (needs SPD matrix)

% eigenvalues and eigenvectors
[V, D] = eig(A);        % A*V = V*D
eigenvalues = diag(D);

% singular value decomposition
[U, S, V] = svd(A);     % A = U*S*V'

% create sparse matrix (memory efficient for large, mostly-zero)
S = sparse(10000, 10000);
S(1, 1) = 5;
S(2, 3) = 10;

% from triplets (row, col, value)
i = [1 2 2 3];
j = [1 1 3 3];
v = [4 5 7 9];
S = sparse(i, j, v, 3, 3);

% convert sparse <-> full
F = full(S);
S2 = sparse(F);

% sparse identity and diagonal
I = speye(1000);
D = spdiags(ones(1000,1), 0, 1000, 1000);

% operations preserve sparsity
nnz(S)                  % number of nonzeros
spy(S)                  % visualize sparsity pattern

A = [1 2 3; 4 5 6];

% reshape (column-major order!)
B = reshape(A, 3, 2);   % [1 5; 4 3; 2 6]
C = reshape(A, 6, 1);   % column vector

% flip and rotate
flipud(A)               % flip up-down
fliplr(A)               % flip left-right
rot90(A)                % rotate 90 degrees
rot90(A, 2)             % rotate 180 degrees

% repmat and repelem
repmat([1 2], 2, 3)     % tile 2x3
repelem([1 2 3], 2)     % [1 1 2 2 3 3]

% concatenation
[A; A]                  % vertical (must have same # cols)
[A A]                   % horizontal (must have same # rows)
horzcat(A, A); vertcat(A, A);

% matrix functions (apply element-wise via fun)
arrayfun(@(x) x^2, A)

% MATLAB implicitly expands (broadcasts) since R2016b
A = [1 2 3; 4 5 6];     % 2x3
b = [10 20 30];         % 1x3
C = A + b;              % 2x3: b added to each row

c = [1; 2];             % 2x1
D = A + c;              % 2x3: c added to each column

% vectorized operations (much faster than loops)
x = linspace(0, 2*pi, 1e6);
y = sin(x) .* cos(x);   % element-wise, no loop needed

% avoid growing arrays in loops
% BAD: for k=1:1e6, s(k) = k; end
% GOOD: preallocate
s = zeros(1, 1e6);
for k = 1:1e6
    s(k) = k;
end

% even better: vectorize
s = 1:1e6;

% bsxfun (legacy, before auto-broadcasting)
E = bsxfun(@plus, A, b);

% least-squares fit (overdetermined system)
A = [1 1; 1 2; 1 3; 1 4];
b = [6; 5; 7; 10];
x = A \ b;             % least-squares solution

% polynomial curve fitting
t = 0:0.1:5;
y = 2*t.^2 + 1 + 0.5*randn(size(t));
p = polyfit(t, y, 2);   % degree-2 fit
yfit = polyval(p, t);

% PCA via SVD
X = randn(100, 5);
Xc = X - mean(X);       % center
[U, S, V] = svd(Xc, 0);
scores = Xc * V;        % principal components
variances = diag(S).^2 / (size(X,1)-1);

% condition number (numerical stability)
cond(A)                 % large = ill-conditioned

x = linspace(0, 2*pi, 200);
y1 = sin(x); y2 = cos(x);

figure;
plot(x, y1, 'b-', 'LineWidth', 2); hold on;
plot(x, y2, 'r--', 'LineWidth', 1.5);
hold off;

% annotations
xlabel('Angle (rad)');
ylabel('Amplitude');
title('Sine and Cosine');
legend('sin(x)', 'cos(x)', 'Location', 'best');
grid on;

% axis control
axis([0 2*pi -1.2 1.2]);
set(gca, 'FontSize', 12, 'FontName', 'Arial');

% multiple subplots
subplot(2, 1, 1); plot(x, y1);
subplot(2, 1, 2); plot(x, y2);

% save figure
saveas(gcf, 'plot.png');
print('-dpng', '-r300', 'plot_hi.png');

[X, Y] = meshgrid(-2:0.1:2);
Z = X .* exp(-X.^2 - Y.^2);

figure;
subplot(2,2,1); surf(X, Y, Z);    title('surf');
subplot(2,2,2); mesh(X, Y, Z);    title('mesh');
subplot(2,2,3); contour(X, Y, Z, 20); title('contour');
subplot(2,2,4); surfc(X, Y, Z);   title('surfc');

% color mapping
colormap(jet(256)); colorbar;
shading interp;        % smooth colors

% 3D line plot (parametric)
t = 0:0.01:10;
plot3(sin(t), cos(t), t, 'LineWidth', 2);
xlabel('x'); ylabel('y'); zlabel('z');

% view angle
view(45, 30);          % azimuth, elevation
axis equal;

data = randn(1000, 1);

% histogram (modern)
histogram(data, 30, 'Normalization', 'pdf');
hold on; x = -4:0.1:4; plot(x, normpdf(x), 'r-', 'LineWidth', 2);

% box plot by group
g = randi([1 3], 100, 1);
boxplot(data, g);

% scatter with marginal histograms
x = randn(200,1); y = x + 0.5*randn(200,1);
scatterhist(x, y);

% bar chart
sales = [10 25 30; 20 15 40];
bar(sales); legend('Q1','Q2','Q3');

% pie chart
pie([30 20 50], {'A','B','C'});

% error bars
y = [1 2 3 4]; e = [0.1 0.2 0.15 0.3];
errorbar(y, e, 'o-');

% animated line plot
figure;
h = plot(NaN, NaN);
axis([0 10 -1 1]); xlabel('t'); ylabel('y');
for t = 0:0.1:10
    set(h, 'XData', [get(h,'XData') t], ...
           'YData', [get(h,'YData') sin(t)]);
    drawnow;
end

% getframe to capture animation
frames = [];
for k = 1:50
    plot(sin(linspace(0, k*pi/10, 200)));
    axis([0 200 -1 1]);
    frames = [frames, getframe];
end
movie(frames, 1, 30);   % play once at 30 fps

% write to video file
v = VideoWriter('anim.mp4', 'MPEG-4');
v.FrameRate = 30; open(v);
for k = 1:50
    plot(sin(linspace(0, k*pi/10, 200)));
    writeVideo(v, getframe);
end
close(v);

% everything in MATLAB graphics is a handle object
h = plot(1:10, 'r-');
set(h, 'LineWidth', 3, 'Marker', 'o', 'MarkerSize', 8);
get(h, 'Color')         % returns [1 0 0]

% find objects
ax = gca;               % current axes
fig = gcf;              % current figure
lines = findobj(ax, 'Type', 'line');

% property exploration
get(ax)                 % list all axes properties
set(ax)                 % list possible values

% text annotations with LaTeX
text(2, 0.5, '\alpha^2 + \beta^2', 'Interpreter', 'tex');
title('$\int_0^1 x^2 dx$', 'Interpreter', 'latex');

% legend with custom location and orientation
legend({'data1','data2'}, 'Location','northwest','Orientation','horizontal');

% export with tight bounding box
exportgraphics(gcf, 'fig.pdf', 'ContentType','vector');

% create a new Simulink model programmatically
mdl = 'my_model';
new_system(mdl);
open_system(mdl);

% add blocks
add_block('simulink/Sources/Sine Wave', [mdl '/Sine']);
add_block('simulink/Sinks/Scope', [mdl '/Scope']);
add_block('simulink/Math Operations/Gain', [mdl '/Gain']);

% set block parameters
set_param([mdl '/Gain'], 'Gain', '2');

% connect blocks (port handles)
sine_out = get_param([mdl '/Sine'], 'PortHandles');
gain_in = get_param([mdl '/Gain'], 'PortHandles');
gain_out = get_param([mdl '/Gain'], 'PortHandles');
scope_in = get_param([mdl '/Scope'], 'PortHandles');
add_line(mdl, sine_out.Outport(1), gain_in.Inport(1));
add_line(mdl, gain_out.Outport(1), scope_in.Inport(1));

% set simulation time and run
set_param(mdl, 'StopTime', '10');
sim(mdl);

% configure solver
set_param(mdl, 'Solver', 'ode45');
set_param(mdl, 'StopTime', '20');
set_param(mdl, 'MaxStep', '0.1');
set_param(mdl, 'RelTol', '1e-4');

% variable-step vs fixed-step
% variable-step: ode45, ode15s, ode23tb (auto-adjusts)
% fixed-step: ode1, ode4, ode5 (real-time/hardware)

% simulate and capture output
out = sim(mdl, 'ReturnWorkspaceOutputs', 'on');
t = out.tout;
y = out.yout;

% simulate from workspace with input
set_param([mdl '/In1'], 'VariableName', 'u');
simIn = Simulink.SimulationInput(mdl);
simIn = setExternalInput(simIn, [t, u]);
out = sim(simIn);

% rapid accelerator for repeated runs
rtp = Simulink.BlockDiagram.buildRapidAcceleratorTarget(mdl);

% create a subsystem
add_block('simulink/Ports & Subsystems/Subsystem', [mdl '/Controller']);

% mask a subsystem to expose parameters
Simulink.Mask.create([mdl '/Controller']);
mask = Simulink.Mask.get([mdl '/Controller']);
mask.addParameter('Type', 'edit', 'Name', 'Kp', 'Prompt', 'Proportional gain');
mask.addParameter('Type', 'edit', 'Name', 'Ki', 'Prompt', 'Integral gain');

% use the mask parameter inside the subsystem
set_param([mdl '/Controller/Gain'], 'Gain', 'Kp');

% library blocks (reusable)
% save model as .slx library: save_system(mdl, 'my_lib.slx');
% set_param(mdl, 'BlockDiagramType', 'Library');

% model reference (separate .slx file)
set_param(mdl, 'ModelReferenceSimulationMode', 'normal');

% variant subsystems (compile-time selection)
add_block('simulink/Ports & Subsystems/Variant Subsystem', [mdl '/Var']);

% Stateflow charts model event-driven logic
% add a Stateflow chart to a model
sf = sfnew;
chart = sfroot.find('-isa', 'Stateflow.Chart');

% MATLAB API to add states and transitions
s1 = Stateflow.State(chart);
s1.Name = 'Off';
s1.Position = [50 50 80 60];
s1.LabelString = 'entry: led = 0;';

s2 = Stateflow.State(chart);
s2.Name = 'On';
s2.Position = [200 50 80 60];
s2.LabelString = 'entry: led = 1;';

% transition from Off to On on event 'press'
t = Stateflow.Transition(chart);
t.Source = s1; t.Destination = s2;
t.Event = 'press';

% simulate chart from MATLAB
% chart inputs: events and data defined in the chart

% generate C/C++ code from a Simulink model
set_param(mdl, 'SystemTargetFile', 'grt.tlc');  % generic real-time
slbuild(mdl);                                     % build

% embed code generation
set_param(mdl, 'SystemTargetFile', 'ert.tlc');   % embedded coder
slbuild(mdl);

% configure for specific hardware
set_param(mdl, 'HardwareImplementation', ...
          'Manufacturer', 'ARM Compatible', 'Type', 'Cortex-M');

% fixed-point conversion (for integer-only hardware)
fxp_cfg = fxptool(mdl);

% generate PLC code (industrial controllers)
plcgeneratecode(mdl);

% deploy to Speedgoat, dSPACE, or other real-time targets
% via vendor-specific blocks and target files

% HIL (hardware-in-the-loop) testing
slrealtime('Speedgoat');

% App Designer (appdesigner command) creates .mlapp files
% Modern replacement for GUIDE (deprecated in R2020a)

% programmatic UI (alternative to App Designer)
f = uifigure('Name', 'My App', 'Position', [100 100 600 400]);

% add components
btn = uibutton(f, 'Text', 'Click Me', ...
    'Position', [50 300 100 30]);
lbl = uilabel(f, 'Text', 'Hello', ...
    'Position', [50 250 200 30]);
fld = uieditfield(f, 'text', 'Position', [50 200 200 30]);

% slider and drop-down
sld = uislider(f, 'Position', [50 150 200 3]);
dd = uidropdown(f, 'Items', {'A','B','C'}, ...
    'Position', [50 100 100 30]);

% axes for plotting
ax = uiaxes(f, 'Position', [300 100 250 250]);
plot(ax, 1:10, rand(1,10));

% callback
btn.ButtonPushedFcn = @(src,event) clickHandler(lbl);

function clickHandler(lbl)
    lbl.Text = ['Clicked at ' datestr(now)];
end

f = uifigure;

% grid layout (responsive)
g = uigridlayout(f, [3 2]);
g.RowHeight = {'1x', '2x', '1x'};
g.ColumnWidth = {'1x', '2x'};

btn1 = uibutton(g); btn1.Layout.Row = 1; btn1.Layout.Column = 1;
btn2 = uibutton(g); btn2.Layout.Row = 1; btn2.Layout.Column = 2;

% tab group
tg = uitabgroup(f);
t1 = uitab(tg, 'Title', 'Plot');
t2 = uitab(tg, 'Title', 'Data');
ax = uiaxes(t1); plot(ax, rand(5));

% panel for grouping
p = uipanel(f, 'Title', 'Settings', 'Position', [20 20 200 150]);

% scrollable container for long content
sf = uiscrollbox(f);
% add many components inside sf

% dialog windows
uigetfile('*.mat', 'Select MAT file');
inputdlg({'Name','Age'}, 'Enter data');
msgbox('Operation complete');

% value-changing callback (slider)
sld = uislider(f);
sld.ValueChangingFcn = @(src,event) onSlider(src, event);

function onSlider(src, event)
    disp(['Dragging: ' num2str(event.Value)]);
end

% value-changed callback (final value)
sld.ValueChangedFcn = @(src,event) disp(['Final: ' num2str(event.Value)]);

% button callback with event data
btn.ButtonPushedFcn = @onButton;
function onButton(src, event)
    fprintf('Button %s pushed
', src.Text);
end

% keyboard callback on figure
f.WindowKeyPressFcn = @(src,event) onKey(event);
function onKey(event)
    fprintf('Key: %s
', event.Key);
end

% listener pattern (any property change)
lh = addlistener(btn, 'Text', 'PostSet', @(src,event) disp('Text changed'));
delete(lh);  % remove listener

% timer for periodic updates
t = timer('ExecutionMode','fixedRate','Period',1, ...
    'TimerFcn',@(~,~) disp('tick'));
start(t);

% package an App Designer app as a .mlappinstall file
% using the Application Compiler (deploytool app)

% programmatically create a packaged app
appFile = matlab.apputil.package(...
    'Name', 'MyPlotTool', ...
    'Summary', 'A simple plotting utility', ...
    'Description', 'Plots user-supplied data.', ...
    'MainFile', 'MyPlotTool.mlapp', ...
    'FolderPath', 'my_app_package');

% install the app
matlab.apputil.install('MyPlotTool.mlappinstall');

% list installed apps
apps = matlab.apputil.getInstalledApps;

% uninstall
matlab.apputil.uninstall('MyPlotTool');

% standalone desktop app (MATLAB Runtime required)
% deploytool -> Application Compiler -> produce .exe

% web app (runs in browser, MATLAB Web App Server)
% deploytool -> Web App

% live plot update from worker
f = uifigure; ax = uiaxes(f); h = plot(ax, NaN, NaN);
for k = 1:100
    set(h, 'XData', 1:k, 'YData', rand(1,k));
    drawnow limitrate;  % throttle for performance
end

% table display
tdata = table({'A';'B';'C'}, [1;2;3], 'VariableNames', {'Name','Val'});
uit = uitable(f, 'Data', tdata.Variables, ...
    'ColumnName', tdata.Properties.VariableNames, ...
    'Position', [20 20 200 150]);

% progress dialog
d = uiprogressdlg(f, 'Title', 'Working', 'Message', 'Processing...');
for k = 1:100
    pause(0.01);
    d.Value = k;
    d.Message = sprintf('%d%% complete', k);
end
close(d);

% context menu
cm = uicontextmenu(f);
m1 = uimenu(cm, 'Text', 'Reset', 'MenuSelectedFcn', @reset);
ax.UIContextMenu = cm;

% read and display an image
img = imread('peppers.png');
imshow(img);
title('RGB image');

% image types
% uint8 (0-255), uint16, double (0-1), logical (binary)
whos img
class(img)              % uint8
size(img)               % rows x cols x channels

% convert between types
d = im2double(img);     % scale to [0,1]
u = im2uint8(d);        % scale to [0,255]
g = im2gray(img);       % RGB to grayscale (R2020b+, replaces rgb2gray)

% write to file
imwrite(g, 'gray.png');
imwrite(g, 'gray.jpg', 'Quality', 90);

% info without loading
info = imfinfo('peppers.png');
disp(info.Width);

% montage of multiple images
fileList = dir('*.png');
montage({fileList.name});

img = im2double(imread('cameraman.tif'));

% Gaussian blur
h = imgaussfilt(img, 2);   % sigma = 2
h2 = imgaussfilt(img, [3 3], 1.5);

% median filter (salt-and-pepper noise)
noisy = imnoise(img, 'salt & pepper', 0.05);
m = medfilt2(noisy, [3 3]);

% custom 2D filter
kernel = fspecial('sobel');
edge_img = imfilter(img, kernel);

% unsharp masking (sharpen)
sharp = imsharpen(img, 'Amount', 0.8, 'Radius', 1.5);

% histogram equalization
eq = histeq(img);          % global
eq2 = adapthisteq(img);    % CLAHE (local)

% image arithmetic
bright = imadd(img, 0.1);
contrast = imadjust(img, [0.3 0.7], []);

bw = imread('text.png');
bw = imbinarize(bw);

% morphological operations
se = strel('disk', 3);
er = imerode(bw, se);       % shrink objects
di = imdilate(bw, se);      % grow objects
op = imopen(bw, se);        % erode then dilate (removes small noise)
cl = imclose(bw, se);       % dilate then erode (fills small holes)

% connected components
cc = bwconncomp(bw);
stats = regionprops(cc, 'Area', 'Centroid', 'BoundingBox');
[areas, idx] = sort([stats.Area], 'descend');

% watershed segmentation
grad = imgradient(bw);
D = -bwdist(~bw);
L = watershed(D);

% edge detection
edges = edge(bw, 'canny');
edges2 = edge(bw, 'sobel');

% thresholding
T = graythresh(img);        % Otsu's method
bw2 = imbinarize(img, T);

img = imread('cameraman.tif');

% corner detection (Harris)
corners = detectHarrisFeatures(img);
[~, strength] = corners;
imshow(img); hold on; plot(corners);

% SURF features (scale-invariant)
points = detectSURFFeatures(img);
[f1, vpts1] = extractFeatures(img, points);

% ORB features (faster, binary)
points = detectORBFeatures(img);

% Hough transform for line detection
bw = edge(img, 'canny');
[H, T, R] = hough(bw);
P = houghpeaks(H, 5);
lines = houghlines(bw, T, R, P);

% image registration
% (find transform between two images)
moving = imread('rotated.png');
fixed = img;
[optimizer, metric] = imregconfig('monomodal');
tform = imregtform(moving, fixed, 'affine', optimizer, metric);
registered = imwarp(moving, tform);

% process many images in a folder
fileList = dir('*.tif');
parfor k = 1:numel(fileList)
    img = imread(fileList(k).name);
    result = processImage(img);
    imwrite(result, ['out_' fileList(k).name]);
end

% block processing for large images
fun = @(block) std2(block);
B = blockproc('large.tif', [256 256], fun);
B = blockproc('large.tif', [256 256], fun, 'Destination', 'out.tif');

% work with very tall/wide images without loading all
% (bigTIFF, multiresolution)
t = Tiff('big.tif', 'r');
imgInfo = t.getImageInfo;

% image datastore for huge collections
imds = imageDatastore('images/', 'IncludeSubfolders', true);
imds.ReadSize = 10;  % read 10 at a time
while hasdata(imds)
    [imgs, info] = read(imds);
    % process batch
end

% tall arrays (out-of-core) for pixel statistics
t = tall(imds);

% unconstrained minimization
f = @(x) (x(1)-2)^2 + (x(2)-3)^2;
x0 = [0 0];
[x, fval] = fminunc(f, x0);

% constrained (with fmincon)
A = []; b = []; Aeq = []; beq = [];
lb = [0 0]; ub = [];
nonlcon = @(x) deal([], x(1) + x(2) - 4);  % x1+x2 <= 4
[x, fval] = fmincon(f, x0, A, b, Aeq, beq, lb, ub, nonlcon);

% linear programming
f_lp = [-1 -1];            % minimize -x1 - x2
A_lp = [1 1; -1 1]; b_lp = [4; 2];
[x, fval] = linprog(f_lp, A_lp, b_lp, [], [], [0;0]);

% quadratic programming
H = [1 -1; -1 2]; c = [-2; -6];
[x, fval] = quadprog(H, c, [], [], [], [], [0;0]);

% set optimization options
opts = optimoptions('fmincon', 'Algorithm','interior-point', ...
    'Display','iter', 'MaxIterations', 200);
[x, fval] = fmincon(f, x0, A, b, Aeq, beq, lb, ub, nonlcon, opts);

% GlobalSearch / MultiStart for non-convex problems
f = @(x) x(1)^2 + x(2)^2 + 10*sin(x(1)) + 10*sin(x(2));
opts = optimoptions('fmincon', 'Algorithm','sqp');
problem = createOptimProblem('fmincon','objective',f, ...
    'x0',[0 0],'lb',[-5 -5],'ub',[5 5],'options',opts);
gs = GlobalSearch;
[x, fval] = run(gs, problem);

% genetic algorithm
[x, fval] = ga(f, 2, [], [], [], [], [-5 -5], [5 5]);

% particle swarm
[x, fval] = particleswarm(f, 2, [-5 -5], [5 5]);

% simulated annealing
[x, fval] = simulannealbnd(f, [0 0], [-5 -5], [5 5]);

% multiobjective (Pareto front)
f_multi = @(x) [x(1)^2 + x(2)^2, (x(1)-1)^2 + (x(2)-1)^2];
[x, fval] = gamultiobj(f_multi, 2, [], [], [], [], [-2 -2], [2 2]);
plot(fval(:,1), fval(:,2), 'o');

% fit a custom model to data
x = linspace(0, 5, 50);
y = 2*exp(-0.5*x) + 0.1*randn(size(x));

% define model and fit
model = fittype('a*exp(-b*x)', 'independent', 'x', ...
    'coefficients', {'a','b'});
f = fit(x', y', model, 'StartPoint', [1 1]);
plot(f, x, y);

% predefined models
f2 = fit(x', y', 'exp2');     % two-term exponential
f3 = fit(x', y', 'poly2');    % quadratic polynomial

% nonlinear least squares with lsqcurvefit
fun = @(p, xdata) p(1)*exp(-p(2)*xdata);
[p, resnorm] = lsqcurvefit(fun, [1 1], x, y);

% confidence intervals
ci = confint(f, 0.95);
predint(f, x', 0.95, 'observation', 'functional');

% smoothing spline
sf = fit(x', y', 'smoothingspline');

% mixed-integer linear programming (TSP-like problems)
f = [-3 -2];             % maximize 3x + 2y
A = [1 1; 2 1]; b = [4; 6];
intcon = 1:2;            % both variables integer
[x, fval] = intlinprog(f, intcon, A, b, [], [], [0;0]);

% assignment problem (Hungarian algorithm)
cost = [4 1 3; 2 3 4; 3 4 1];
[assignment, cost_total] = matchpairs(cost, false);

% traveling salesman via intlinprog
n = 8;  % cities
distances = rand(n, n); distances = distances + distances';
% (build subtour elimination constraints iteratively)

% knapsack
weights = [2 3 4 5]; values = [3 4 5 6];
W = 5;
f = -values;
A = weights; b = W;
intcon = 1:4;
[x, fval] = intlinprog(f, intcon, A, b, [], [], zeros(4,1), ones(4,1));

% surrogate optimization (expensive black-box)
obj = @(x) expensive_eval(x);
[x, fval] = surrogateopt(obj, lb, ub);

% 1. choose solver based on problem type
%    smooth? -> fminunc/fmincon
%    linear? -> linprog/intlinprog
%    integer? -> intlinprog
%    nonsmooth/global? -> ga/particleswarm/GlobalSearch

% 2. provide gradients and Jacobians when possible
f = @(x) x(1)^2 + x(2)^2;
gradf = @(x) [2*x(1); 2*x(2)];
prob = optimproblem('Objective', f);
prob.Objective.Gradient = gradf;

% 3. scale variables to ~[1, 1]
x_scaled = (x - x0) / scale;

% 4. check optimality and feasibility
output.exitflag       % 1 = converged
output.firstorderopt  % should be near 0

% 5. warm start sequential problems
opts = optimoptions('quadprog','ObjectiveLimit',1e-8);
x_new = quadprog(H, c, A, b, [], [], lb, ub, x_old, opts);

% 6. use problem-based framework (cleaner)
x = optimvar('x', 2, 'LowerBound', 0);
prob = optimproblem('Objective', x(1)^2 + x(2)^2);
prob.Constraints.c1 = x(1) + x(2) >= 1;
sol = solve(prob);

% start a parallel pool (workers = CPU cores)
pool = gcp('nocreate');
if isempty(pool)
    parpool;  % start default pool
end

% parfor: parallel for loop
n = 1e6;
results = zeros(1, n);
parfor k = 1:n
    results(k) = expensive_computation(k);
end

% restrictions:
% - iterations must be independent (no data dependencies)
% - variables classified as loop, sliced, broadcast, reduction
% - cannot break/continue

% reduction variables (combine results)
total = 0;
parfor k = 1:n
    total = total + compute(k);
end

% sliced variables (each iteration writes distinct element)
out = zeros(n, 1);
parfor k = 1:n
    out(k) = k^2;
end

% control pool size
delete(gcp('nocreate'));  % shut down
parpool('local', 4);      % 4 workers

% spmd: single program, multiple data
spmd
    labidx              % worker index (1 to n)
    numlabs             % total workers
    A = rand(1000) * labidx;
    % each lab has its own A
end

% access Composite (per-worker data)
A{1}                    % A from worker 1
A{2}                    % A from worker 2

% distributed arrays (spread across workers)
D = distributed.rand(10000);
size(D)                 % 10000x10000 total
L = localpart(D);       % this worker's chunk

% communicate between labs
spmd
    data = rand(100, 1);
    % send to next lab, receive from previous
    next_lab = mod(labidx, numlabs) + 1;
    prev_lab = mod(labidx - 2, numlabs) + 1;
    rcvd = labSendrecv(data, next_lab, prev_lab);
end

% gather distributed to local
D = distributed.rand(1000);
full_matrix = gather(D);

% move data to GPU
A = gpuArray(rand(10000));
B = gpuArray(eye(10000));

% operations execute on GPU
C = A * B;              % GPU matrix multiply
D = A + B;
E = sum(A, 1);

% gather results back to CPU
C_cpu = gather(C);

% many built-in functions are GPU-enabled
x = gpuArray.linspace(-10, 10, 1e6);
y = sin(x) .* exp(-x.^2);

% arrayfun on GPU (custom element-wise functions)
f = @(a, b) a.^2 + b.^2;
result = arrayfun(f, A, B);

% check GPU
gpuDevice            % info about current GPU
gpuDeviceCount       % number of available GPUs

% page multiple GPUs
parfor k = 1:4
    gpuDevice(k);
    A = gpuArray(rand(5000));
    % ... compute on GPU k
end

% run a function in background
job = batch(@my_function, 0, {arg1, arg2});

% check status
job.State              % 'queued', 'running', 'finished'
wait(job);             % block until done
load(job);             % retrieve results

% batch with parfor inside
job = batch(@my_parfor_func, 1, {input_data}, 'Pool', 4);

% job and task objects (lower-level)
job = createJob(pool);
task = createTask(job, @my_func, 1, {arg1});
submit(job);
wait(job);
result = fetchOutputs(task);

% cluster types
% 'local' - on this machine
% 'MJS' - MATLAB Job Scheduler (server cluster)
% 'SLURM'/'PBS' - HPC scheduler integration
c = parcluster('local');
c.NumWorkers = 8;
c.saveProfile;

% profile parallel code
profile on -parallel
parfor k = 1:100
    compute(k);
end
profile off
profview                % visualize

% measure speedup
tic;
serial_result = serial_compute();
t_serial = toc;

tic;
parfor_result = parallel_compute();
t_parallel = toc;

speedup = t_serial / t_parallel;
fprintf('Speedup: %.2fx using %d workers
', speedup, pool.NumWorkers);

% avoid common pitfalls:
% 1. too little work per iteration (overhead dominates)
% 2. large broadcast variables (transfer cost)
% 3. false sharing via reduction variables

% data queue for live updates from workers
Q = parallel.pool.DataQueue;
afterEach(Q, @update_plot);
parfor k = 1:n
    % ... compute ...
    send(Q, partial_result);
end

% class definition in a file named Point.m
classdef Point
    properties
        x = 0       % public, with default
        y = 0
    end
    properties (SetAccess = protected)
        id          % settable only by class methods
    end
    properties (Constant)
        PI = 3.14159
    end

    methods
        function obj = Point(x, y)
            if nargin > 0
                obj.x = x;
                obj.y = y;
            end
            obj.id = randi(1e6);
        end

        function r = distance(obj)
            r = sqrt(obj.x^2 + obj.y^2);
        end

        function disp(obj)
            fprintf('Point(%.2f, %.2f) [id=%d]
', ...
                obj.x, obj.y, obj.id);
        end
    end
end

% usage
p = Point(3, 4);
p.distance()           % 5
disp(p)

% value class (default): copies on assignment
classdef Vec
    properties
        data
    end
    methods
        function obj = Vec(d)
            obj.data = d;
        end
    end
end

v1 = Vec([1 2 3]);
v2 = v1;             % COPY
v2.data = [4 5 6];
disp(v1.data)        % still [1 2 3]

% handle class: references shared
classdef HVec < handle
    properties
        data
    end
    methods
        function obj = HVec(d)
            obj.data = d;
        end
    end
end

h1 = HVec([1 2 3]);
h2 = h1;             % REFERENCE (same object)
h2.data = [4 5 6];
disp(h1.data)        % [4 5 6] - changed!

% handle classes support delete method
% isequal(h1, h2) -> true (same object)

% single inheritance
classdef Dog < Animal
    methods
        function obj = Dog(name)
            obj@Animal(name);  % call superclass constructor
        end
        function sound(obj)
            disp([obj.name ' says Woof']);
        end
    end
end

% multiple inheritance
classdef FlyingFish < Fish & Bird
    methods
        function obj = FlyingFish(name)
            obj@Fish(name);
            obj@Bird(name);
        end
    end
end

% abstract methods and classes
classdef Shape
    methods (Abstract)
        area(obj)
        perimeter(obj)
    end
end

classdef Circle < Shape
    properties
        radius
    end
    methods
        function obj = Circle(r)
            obj.radius = r;
        end
        function a = area(obj)
            a = pi * obj.radius^2;
        end
        function p = perimeter(obj)
            p = 2 * pi * obj.radius;
        end
    end
end

% polymorphism
shapes = {Circle(2), Square(3)};
for k = 1:numel(shapes)
    fprintf('Area: %.2f
', shapes{k}.area());
end

classdef TemperatureSensor < handle
    properties
        Temperature = 20
    end
    events
        TemperatureChanged
        Overheated
    end
    methods
        function set.Temperature(obj, t)
            obj.Temperature = t;
            notify(obj, 'TemperatureChanged');
            if t > 100
                notify(obj, 'Overheated');
            end
        end
        function obj = TemperatureSensor()
            % nothing
        end
    end
end

% create sensor and listen
s = TemperatureSensor();
addlistener(s, 'TemperatureChanged', @onTempChanged);
addlistener(s, 'Overheated', @onOverheat);

function onTempChanged(src, event)
    fprintf('Temp now: %.1f
', src.Temperature);
end

function onOverheat(src, event)
    warning('Overheated at %.1f!', src.Temperature);
end

% trigger
s.Temperature = 25;     % prints "Temp now: 25.0"
s.Temperature = 105;    % prints temp + warning

classdef Complex
    properties
        re
        im
    end
    methods
        function obj = Complex(re, im)
            obj.re = re; obj.im = im;
        end
        function r = plus(a, b)            % a + b
            r = Complex(a.re + b.re, a.im + b.im);
        end
        function r = mtimes(a, b)          % a * b
            r = Complex(a.re*b.re - a.im*b.im, ...
                        a.re*b.im + a.im*b.re);
        end
        function r = abs(obj)              % |obj|
            r = sqrt(obj.re^2 + obj.im^2);
        end
        function disp(obj)
            fprintf('%.2f %+.2fi
', obj.re, obj.im);
        end
    end
    methods
        function r = subsref(obj, S)       % obj(...)
            if strcmp(S(1).type, '()')
                r = obj.re + 1i*obj.im;
            end
        end
        function obj = subsasgn(obj, S, val)  % obj(...) = val
            % custom assignment
        end
    end
end

c1 = Complex(1, 2); c2 = Complex(3, 4);
c3 = c1 + c2;          % uses plus
c4 = c1 * c2;          % uses mtimes
abs(c1)                % uses abs

% set breakpoints in the Editor or programmatically
dbstop in my_function at 42      % line 42 of my_function.m
dbstop in my_function at 42 if x > 10  % conditional
dbstop if error                  % pause on any error
dbstop if naninf                 % pause on NaN/Inf
dbstop if warning                % pause on warnings
dbclear all                      % clear all breakpoints
dbclear in my_function           % clear in one file

% at a breakpoint, the command window enters debug mode K>>
whos                             % inspect workspace
x                                % view variables
dbstep                           % step to next line
dbstep in                        % step into function
dbstep out                       % step out of function
dbcont                           % continue execution
dbstack                          % show call stack
dbup / dbdown                    % move up/down stack frames
dbquit                           % exit debug mode

% error recovery
try
    risky_operation();
catch ME
    fprintf('Error: %s
', ME.message);
    fprintf('Stack:
');
    for k = 1:numel(ME.stack)
        fprintf('  %s (line %d)
', ME.stack(k).name, ME.stack(k).line);
    end
end

% profile a script or function
profile on
my_expensive_script();
profile off

% view results
profile viewer                 % GUI
p = profile('info');           % structure
profsave(p, 'profile_results'); % HTML report

% time individual operations
tic; A = rand(5000); A = A * A; toc;   % ~0.5s

% timeit for accurate single-function timing
f = @() sort(rand(1e6,1));
t = timeit(f);

% memory profiling
profile on -memory
% ... run code ...
profile viewer   % shows memory allocation

% identify bottlenecks
% 1. profile viewer shows time per function/line
% 2. focus on the top few hottest lines
% 3. vectorize, preallocate, or use MEX for those

% BAD: growing array in loop
tic;
s = [];
for k = 1:100000
    s = [s, k^2];     % reallocates every iteration!
end
toc;                  % ~10 seconds

% GOOD: preallocate
tic;
s = zeros(1, 100000);
for k = 1:100000
    s(k) = k^2;
end
toc;                  % ~0.01 seconds

% BEST: vectorize
tic;
s = (1:100000).^2;
toc;                  % ~0.001 seconds

% vectorize conditional logic
x = rand(10000, 1);
% BAD:
% for k = 1:numel(x)
%     if x(k) > 0.5, y(k) = 1; else, y(k) = 0; end
% end
% GOOD:
y = double(x > 0.5);

% vectorize with logical indexing
y = zeros(size(x));
y(x > 0.5) = 1;
y(x > 0.8) = 2;

% use arrayfun/gpuArray only when truly needed
% (loops with preallocation are often faster than arrayfun)

% call C/C++/Fortran code from MATLAB via MEX
% write a C file (my_func.c):
/*
#include "mex.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
    double *x = mxGetPr(prhs[0]);
    size_t n = mxGetNumberOfElements(prhs[0]);
    plhs[0] = mxCreateDoubleMatrix(n, 1, mxREAL);
    double *y = mxGetPr(plhs[0]);
    for (size_t k = 0; k < n; k++) y[k] = x[k] * x[k];
}
*/

% compile
mex my_func.c           % produces my_func.mexw64
% usage:
y = my_func([1 2 3]);   % [1 4 9]

% C++ with MATLAB Data API (modern, R2018a+)
% #include "mex.hpp", "mexAdapter.hpp"
% use matlab::data::TypedArray<double>

% call MATLAB from C/C++ application (MATLAB Engine API)
% #include "engine.h"
% eng = engOpen("");
% engPutVariable(eng, "x", x);
% engEvalString(eng, "y = sqrt(x)");

% loadlibrary for generic shared libraries
loadlibrary('mylib.dll', 'mylib.h');
calllib('mylib', 'my_function', arg);
unloadlibrary('mylib');

% check memory usage
[usr, sys] = memory;
disp(usr.PhysicalMemory.Available);

% view variable sizes
whos                  % all variables with sizes
[x_size, x_bytes] = compute_size(x);

% clear large variables when done
clear big_array;

% use single precision when possible
A = rand(5000, 'single');   % half the memory of double
B = A * A;                   % single precision matmul

% sparse for mostly-zero matrices
S = sparse(1e5, 1e5);
S(1, 1) = 1;

% memory-mapped files for huge data
m = memmapfile('large.bin', 'Format', 'double', 'Offset', 0, ...
    'Repeat', 1e8);
m.Data(1)              % access without loading all

% tall arrays for out-of-core
ds = datastore('big.csv');
t = tall(ds);
mean_t = mean(t.Var1);
gather(mean_t);        % triggers computation

% avoid copies: in-place operations
A = A + 1;             % may modify A in place (since R2018b)
A(:) = A(:) + 1;       % definitely in place

% tiledlayout (R2019b+) replaces subplot
tiledlayout(2, 3, 'TileSpacing', 'compact', 'Padding', 'compact');

nexttile; plot(1:10);
nexttile; scatter(rand(1, 50), rand(1, 50));
nexttile; bar(rand(5, 3));
nexttile; imagesc(magic(10));
nexttile; histogram(randn(1000, 1), 30);
nexttile; pie([30 20 50]);

% span multiple tiles
nexttile(7, [1 2]); plot(sin(0:0.1:2*pi));

% shared title and labels
title('Dashboard Demo');
xlabel('X'); ylabel('Y');

% legacy subplot
subplot(2, 2, 1); plot(1:10);
subplot(2, 2, 2); plot(10:-1:1);

x = linspace(0, 2*pi, 100);
plot(x, sin(x));

% axis limits
xlim([0 2*pi]); ylim([-1.2 1.2]);

% ticks and labels
xticks(0:pi/2:2*pi);
xticklabels({'0', '\pi/2', '\pi', '3\pi/2', '2\pi'});

% grid and minor grid
grid on; grid minor on;

% axis appearance
box on;
ax = gca;
ax.FontSize = 12;
ax.FontName = 'Courier';
ax.LineWidth = 1.5;
ax.XColor = [0.5 0.5 0.5];
ax.YDir = 'reverse';

% log scale
semilogy(1:100, 2.^(1:100));
set(gca, 'YScale', 'log');

x = linspace(0, 2*pi, 100);
y1 = sin(x); y2 = cos(x); y3 = sin(x) + cos(x);

% hold on for multiple lines
plot(x, y1, 'r-', 'LineWidth', 2); hold on;
plot(x, y2, 'b--', 'LineWidth', 2);
plot(x, y3, 'g:', 'LineWidth', 2); hold off;

% legend with location
legend('sin(x)', 'cos(x)', 'sin+cos', ...
       'Location', 'best', 'Orientation', 'horizontal');

% legend with TeX
legend('$\sin(x)$', '$\cos(x)$', 'Interpreter', 'latex');

% specify which lines to include
h1 = plot(x, y1); hold on;
h2 = plot(x, y2);
h3 = plot(x, y3);
legend([h1 h3], 'sin', 'sin+cos');

[X, Y] = meshgrid(-2:0.1:2, -2:0.1:2);
Z = X .* exp(-X.^2 - Y.^2);

% surface plot
figure;
surf(X, Y, Z);
colorbar; shading interp;

% mesh plot (wireframe)
figure;
mesh(X, Y, Z);
hidden off;                      % show hidden lines

% contour
figure;
contour(X, Y, Z, 20); colorbar;
contourf(X, Y, Z, 20); colorbar; % filled
[c, h] = contour(X, Y, Z, 20); clabel(c, h);

% 3D parametric
t = 0:0.01:10;
plot3(sin(t), cos(t), t);
xlabel('x'); ylabel('y'); zlabel('z');
view(45, 30);                    % azimuth, elevation

% animate a sine wave
figure;
h = plot(NaN, NaN);
xlim([0 2*pi]); ylim([-1 1]);

for k = 1:100
  x = linspace(0, 2*pi, k);
  set(h, 'XData', x, 'YData', sin(x));
  drawnow;
end

% capture frames and save as video
v = VideoWriter('animation.mp4', 'MPEG-4');
v.FrameRate = 30; open(v);

for k = 1:60
  plot(sin(0:0.1:k/10));
  frame = getframe(gcf);
  writeVideo(v, frame);
end
close(v);

% getframe for static snapshots
figure; plot(1:10);
frame = getframe;                % struct with cdata

A = [1 2; 3 4];
b = [5; 6];

% solve Ax = b
x = A \ b;                       % preferred (LU-based)
x = inv(A) * b;                  % avoid — slower, less stable

% factorizations
[L, U, P] = lu(A);               % PA = LU
[Q, R] = qr(A);                  % A = QR
[V, D] = eig(A);                 % A V = V D
[U, S, V] = svd(A);              % A = U S V'

% properties
det(A); trace(A); rank(A); cond(A);
null(A);                         % null space
orth(A);                         % orthonormal basis of range

% matrix functions
expm(A); logm(A); sqrtm(A);      % matrix (not element-wise)

A = [1 2; 3 4];
B = [5 6; 7 8];

% element-wise (use .)
A .* B                           % [5 12; 21 32]
A ./ B                           % element-wise divide
A .^ 2                           % [1 4; 9 16]
A .^ B                           % element-wise power

% matrix ops (no .)
A * B                            % matrix product
A / B                            % A * inv(B)
A \ B                            % inv(A) * B
A^2                              % A * A (matrix power)
A'                               % conjugate transpose
A.'                              % non-conjugate transpose

% common gotcha
x = [1 2 3];
x * x                            % ERROR: inner dims mismatch
x * x'                           % 14 (dot product)
x .* x                           % [1 4 9] (element-wise)

A = magic(4);                  % 4x4 matrix

% linear indexing (column-major)
A(5)                             % A(1, 2)
A(:)                             % column vector of all elements
A(2:5)                           % elements 2-5

% submatrix
A(1:2, 2:4)                      % rows 1-2, cols 2-4
A(:, 3)                          % all rows, col 3
A(end, :)                        % last row
A(1:2:end, :)                    % odd rows (stride)

% logical indexing
mask = A > 10;
A(mask)                          % flat list of values > 10
A(A > 10) = 0;                   % threshold in-place

% reshape, permute, flip
B = reshape(A, 2, 8);
C = permute(A, [2 1]);           % transpose N-D
flip(A); flipud(A); fliplr(A);
repmat(A, 2, 3);                 % tile 2x3

% create sparse
S = sparse(1000, 1000);
S(1, 1) = 5;
S(500, 500) = 10;

% from triplets
i = [1 2 3 4];
j = [1 2 3 4];
v = [10 20 30 40];
S = sparse(i, j, v, 4, 4);

% sparse identity and diagonal
S = speye(1000);
S = spdiags(ones(1000, 1), 0, 1000, 1000);

% operations stay sparse
S = S + S';
S = S * S;
x = S \ b;                       % uses sparse solver

% convert
F = full(S);                     % sparse -> full
S = sparse(F);                   % full -> sparse
nnz(S)                           % count non-zeros
spy(S)                           % visualize sparsity

% R2016b+ — automatic broadcasting
A = magic(3);                    % 3x3
b = [1 2 3];                     % 1x3

% older MATLAB: error (dims mismatch)
% modern MATLAB: b is broadcast across rows
C = A + b;                       % 3x3 result
C = A - b';
C = A .* b;

% common uses
x = (1:5)';                      % column
y = (1:5);                       % row
outer = x * y;                   % 5x5 outer product
grid = x + y;                    % 5x5 addition grid

% mean across rows
M = magic(5);
row_means = mean(M, 2);          % 5x1
centered = M - row_means;        % broadcast subtraction

% bsxfun (legacy, pre-R2016b)
C = bsxfun(@plus, A, b);         % same as A + b

% heterogeneous container
C = {'hello', 42, [1 2 3], magic(3)};

% access
C{1}                             % 'hello' (contents)
C(1)                             % {1x1 cell} (cell)
C{4}(2, 2)                       % access inside

% build dynamically
C = {};
for k = 1:5
  C{k} = rand(k);
end

% multi-element
[a, b, c] = C{1:3};              % unpack
C{1:3}                           % comma-separated list

% cellfun
nums = {1, 2, 3, 4};
squared = cellfun(@(x) x^2, nums);
lengths = cellfun(@length, C);

% convert
cell([1 2 3])                    % {1, 2, 3}
cell2mat({1 2; 3 4})             % [1 2; 3 4]
mat2cell(magic(4), [2 2], [2 2]) % 2x2 cell of blocks

% scalar struct
s.name = 'Alice';
s.age = 30;
s.scores = [90 85 92];

% struct constructor
s = struct('name', 'Bob', 'age', 25, 'scores', {[80 70]});

% struct array
students(1).name = 'Alice';
students(1).age = 30;
students(2).name = 'Bob';
students(2).age = 25;

% access across array
names = {students.name};         % all names (cell)
ages = [students.age];           % all ages (vector)

% nested
s.course.code = 'CS101';
s.course.credits = 3;
s.course.instructor.name = 'Dr. Smith';

% field operations
fieldnames(s)                    % list fields
isfield(s, 'age')                % check field
rmfield(s, 'age')                % remove field
orderfields(s)                   % sort fields

% create table
Names = {'Alice'; 'Bob'; 'Carol'};
Age = [30; 25; 42];
Score = [90; 80; 95];
T = table(Names, Age, Score);

% access
T.Age                            % column as vector
T.(2)                            % column by index
T{1, 2}                          % row 1, col 2 (cell-style)
T(1:2, :)                        % first 2 rows

% add/remove columns
T.Grade = {'A'; 'B'; 'A'};
T.Grade = [];                    % remove column

% filtering
T(T.Age > 26, :)                 % rows where Age > 26
find(T.Age > 26);

% sorting
T = sortrows(T, 'Age');
T = sortrows(T, {'Age', 'Score'}, 'descend');

% summary
summary(T)
T.Properties.VariableNames

% read/write
T = readtable('data.csv');
writetable(T, 'out.csv')

% create timetable with timestamps
Time = datetime(2024, 1, 1) + hours(0:23)';
Temp = 20 + 5 * randn(24, 1);
TT = timetable(Time, Temp);

% indexing by time
TT.timerange(datetime(2024,1,1,6,0,0), datetime(2024,1,1,12,0,0))
TT(datetime(2024,1,1,10,0,0), :) % row at specific time
TT('2024-01-01 10:00:00', :)

% resample / synchronize
TT_30min = retime(TT, 'regular', 'mean', 'TimeStep', minutes(30));
TT_hourly = retime(TT, 'hourly', 'mean');

% aggregate
TT_daily = retime(TT, 'daily', @mean);
TT_daily_max = retime(TT, 'daily', 'max');

% combine timetables
TT2 = synchronize(TT1, TT2, 'union');

% lagging
TT.PrevTemp = lag(TT.Temp, 1);

% containers.Map (pre-R2022b)
m = containers.Map;
m('apple') = 1;
m('banana') = 2;
m('cherry') = 3;

% access
m('apple')                       % 1
m.Keys                           % {'apple', 'banana', 'cherry'}
m.Values                         % [1, 2, 3]
isKey(m, 'apple')                % true
remove(m, 'banana');

% iterate
keys = m.keys;
for k = 1:length(keys)
  fprintf('%s -> %d\n', keys{k}, m(keys{k}));
end

% R2022b+ dictionary (faster, typed)
d = dictionary({'apple', 'banana'}, [1, 2]);
d('cherry') = 3;
d('apple')                        % 1
keys(d); values(d)

% minimize objective with constraints
% min f(x) s.t. A*x <= b, Aeq*x = beq, lb <= x <= ub, c(x) <= 0, ceq(x) = 0

fun = @(x) (x(1) - 2)^2 + (x(2) - 3)^2;

% linear inequality: x1 + x2 <= 4
A = [1 1]; b = 4;
% linear equality: x1 - x2 = 0
Aeq = [1 -1]; beq = 0;
% bounds
lb = [0 0]; ub = [5 5];
% nonlinear: x1^2 + x2^2 >= 1 (so -(x1^2 + x2^2 - 1) <= 0)
nonlcon = @(x) deal(-(x(1)^2 + x(2)^2 - 1), []);

x0 = [1 1];
[x, fval, exitflag, output] = fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon);

% options
opts = optimoptions('fmincon', 'Display', 'iter', ...
                    'Algorithm', 'interior-point');
[x, fval] = fmincon(fun, x0, A, b, Aeq, beq, lb, ub, nonlcon, opts);

% unconstrained minimization
fun = @(x) 100*(x(2) - x(1)^2)^2 + (1 - x(1))^2;  % Rosenbrock

x0 = [-1.2 1];
[x, fval, flag, out] = fminunc(fun, x0);

% with gradient (faster, more accurate)
syms x1 x2
f = 100*(x2 - x1^2)^2 + (1 - x1)^2;
grad = gradient(f, [x1, x2]);
fun_grad = matlabFunction(f, grad, 'Vars', {x1, x2});

opts = optimoptions('fminunc', 'SpecifyObjectiveGradient', true, ...
                    'Display', 'iter');
[x, fval] = fminunc(@(x) fun_grad(x(1), x(2)), x0, opts);

% check optimality
disp(out.firstorderopt)           % should be near 0

% Hessian option
opts.Hessian = 'bfgs';            % quasi-Newton (default)
opts.Hessian = 'lbfgs';           % for large sparse problems

% min f'*x s.t. A*x <= b, Aeq*x = beq, lb <= x <= ub

% example: min -x1 - 2*x2
%   s.t.  x1 + x2 <= 4
%          x1 + 3*x2 <= 6
%          x1, x2 >= 0
f = [-1; -2];
A = [1 1; 1 3];
b = [4; 6];
lb = [0; 0];

opts = optimoptions('linprog', 'Display', 'iter');
[x, fval, exitflag, output] = linprog(f, A, b, [], [], lb, [], opts);

% integer programming (intlinprog)
intcon = 1:2;                     % variables 1 and 2 are integers
[x, fval] = intlinprog(f, intcon, A, b, [], [], lb, []);

% binary variables: set lb=0, ub=1, and add to intcon
lb = [0; 0]; ub = [1; 1];
[x, fval] = intlinprog(f, intcon, A, b, [], [], lb, ub);

% fit y = a*exp(b*x) to data
xdata = linspace(0, 5, 50)';
ydata = 2.5 * exp(-0.7 * xdata) + 0.05 * randn(50, 1);

model = @(p, x) p(1) * exp(p(2) * x);

x0 = [1, -1];
[params, resnorm, residual, exitflag] = lsqcurvefit(model, x0, xdata, ydata);

% with bounds
lb = [0, -5]; ub = [10, 0];
[params, resnorm] = lsqcurvefit(model, x0, xdata, ydata, lb, ub);

% polyfit for polynomials
p = polyfit(xdata, ydata, 3);     % cubic
yfit = polyval(p, xdata);

% fit with custom output
opts = optimoptions('lsqcurvefit', 'Display', 'off', ...
                    'Algorithm', 'levenberg-marquardt');
[params, resnorm] = lsqcurvefit(model, x0, xdata, ydata, [], [], opts);

% goodness of fit
sse = sum((ydata - model(params, xdata)).^2);
sst = sum((ydata - mean(ydata)).^2);
r2 = 1 - sse / sst;

% MultiStart: run local solver from many points
fun = @(x) x(1)^2 + x(2)^2 + 10*sin(x(1)) + 10*sin(x(2));
ms = MultiStart('Display', 'iter');
opts = optimoptions('fmincon', 'Algorithm', 'interior-point');
problem = createOptimProblem('fmincon', 'objective', fun, ...
                            'x0', [0 0], 'lb', [-5 -5], 'ub', [5 5]);
[x, fval] = run(ms, problem, 50);  % 50 starting points

% GlobalSearch: smarter, fewer starts
gs = GlobalSearch;
[x, fval] = run(gs, problem);

% ga (genetic algorithm)
[x, fval] = ga(fun, 2, [], [], [], [], [-5 -5], [5 5]);

% simulannealbnd (simulated annealing)
[x, fval] = simulannealbnd(fun, [0 0], [-5 -5], [5 5]);

% patternsearch
[x, fval] = patternsearch(fun, [0 0]);

% read and display
img = imread('photo.jpg');
imshow(img);

% info
size(img)                         % rows x cols x channels
class(img)                        % usually uint8
imfinfo('photo.jpg')              % metadata

% convert types
gray = rgb2gray(img);             % RGB -> gray
img_double = im2double(img);      % uint8 [0,255] -> double [0,1]
img_uint8 = im2uint8(img_double);

% write
imwrite(gray, 'gray.png');
imwrite(img, 'compressed.jpg', 'Quality', 80);

% display multiple
figure;
subplot(1, 2, 1); imshow(img); title('Original');
subplot(1, 2, 2); imshow(gray); title('Grayscale');

% pixel info
impixelinfo;                      % interactive pixel values
img(100, 200, :)                  % pixel at (row, col)

img = imread('cameraman.tif');
img = im2double(img);

% Gaussian blur
blurred = imgaussfilt(img, 2);    % sigma = 2
blurred = imgaussfilt(img, 2, 'FilterSize', 7);

% median filter (great for salt-and-pepper noise)
noisy = imnoise(img, 'salt & pepper', 0.05);
denoised = medfilt2(noisy, [3 3]);

% custom kernel
h = fspecial('sobel');            % edge detection
edges = imfilter(img, h);
h = fspecial('laplacian', 0);
sharpened = img - imfilter(img, h);

% 2D convolution vs correlation
out_conv = imfilter(img, h, 'conv');
out_corr = imfilter(img, h, 'corr');  % default

% unsharp masking
sharp = imsharpen(img, 'Amount', 1.5, 'Radius', 2);

bw = imread('text.png');
bw = imbinarize(bw);              % grayscale -> binary

% basic ops
eroded = imerode(bw, strel('square', 3));
dilated = imdilate(bw, strel('disk', 2));
opened = imopen(bw, strel('square', 5));
closed = imclose(bw, strel('square', 5));

% structuring elements
se = strel('square', 3);
se = strel('disk', 5);
se = strel('line', 10, 45);       % length 10, angle 45 deg
se = strel('octagon', 3);

% remove small objects
clean = bwareaopen(bw, 50);       % remove <50 pixels

% skeleton
skel = bwmorph(bw, 'skel', Inf);

% fill holes
filled = imfill(bw, 'holes');

% boundary
boundary = bwperim(bw);

img = imread('coins.png');
img = im2double(img);

% edges
edges_sobel = edge(img, 'sobel');
edges_canny = edge(img, 'canny', [0.1 0.2], 1.5);
edges_prewitt = edge(img, 'prewitt');

% thresholding
bw = imbinarize(img);
bw_otsu = imbinarize(img, 'adaptive');  % local threshold

% connected components
cc = bwconncomp(bw);
stats = regionprops(cc, 'Area', 'Centroid', 'BoundingBox');
areas = [stats.Area];
[~, idx] = sort(areas, 'descend');
biggest = stats(idx(1));

% watershed segmentation
grad = imgradient(img);
L = watershed(grad);
L(L == 0) = NaN;
imshow(img); hold on; imshow(L, []);

img = imread('peppers.png');

% color space conversion
hsv = rgb2hsv(img);
lab = rgb2lab(img);
gray = rgb2gray(img);
rgb = hsv2rgb(hsv);

% extract channels
R = img(:, :, 1);
G = img(:, :, 2);
B = img(:, :, 3);

% Fourier transform
f = fft2(double(gray));
f = fftshift(f);                  % center zero freq
mag = log(1 + abs(f));
imshow(mag, []);

% inverse
f_filtered = f;
f_filtered(abs(f) < 100) = 0;    % low-pass filter
filtered = real(ifft2(ifftshift(f_filtered)));
imshow(filtered, []);

% resize, rotate
small = imresize(img, 0.5);
rotated = imrotate(img, 45, 'bilinear', 'crop');
cropped = imcrop(img, [50 50 200 150]);  % [x y w h]

% Open App Designer:
%   >> appdesigner
% Drag components from the left, set properties on the right.

% Callback structure (auto-generated):
function ButtonPushed(app, event)
    x = linspace(0, 2*pi, 100);
    plot(app.UIAxes, x, sin(x));
    app.ResultEditField.Value = 'Plot done';
end

% Access components via app.<ComponentName>
%   app.UIAxes, app.Button, app.EditField, app.Slider, ...

% Property panel: add custom properties
properties (Access = private)
    data;                          % shared between callbacks
end

% Startup callback:
function startupFcn(app)
    app.data = rand(10);
    app.Slider.Value = 50;
end

% Save as .mlapp, run with:
%   >> myApp  (if named myApp.mlapp)
% Or export to a standalone .m file.

% modern uifigure (R2016b+) — preferred
fig = uifigure('Name', 'My App', 'Position', [100 100 600 400]);

ax = uiaxes(fig, 'Position', [200 100 350 250]);
plot(ax, 1:10);

btn = uibutton(fig, 'Text', 'Plot', 'Position', [50 100 100 30]);
btn.ButtonPushedFcn = @(~,~) plot(ax, rand(1, 10));

slider = uislider(fig, 'Position', [50 200 100 3]);
slider.ValueChangedFcn = @(s,~) disp(s.Value);

% dropdown
dd = uidropdown(fig, 'Items', {'a', 'b', 'c'}, 'Position', [50 250 100 30]);
dd.ValueChangedFcn = @(s,~) fprintf('Selected: %s\n', s.Value);

% list box
lb = uilistbox(fig, 'Items', {'one', 'two', 'three'}, ...
               'Multiselect', 'on', 'Position', [50 50 100 80]);

% wait for user to close
uiwait(fig);

fig = uifigure;

% label
lbl = uilabel(fig, 'Text', 'Enter value:', 'Position', [20 200 100 22]);

% edit field
ef = uieditfield(fig, 'numeric', 'Value', 42, 'Position', [130 200 100 22]);
ef.ValueChangedFcn = @(s,~) fprintf('New: %g\n', s.Value);

% button group with radio buttons
bg = uibuttongroup(fig, 'Title', 'Mode', 'Position', [20 100 200 80]);
r1 = uiradiobutton(bg, 'Text', 'A', 'Position', [10 40 100 22]);
r2 = uiradiobutton(bg, 'Text', 'B', 'Position', [10 10 100 22]);
r1.Value = true;

% check box
cb = uicheckbox(fig, 'Text', 'Enable', 'Position', [20 50 100 22]);

% table
tdata = table({'A'; 'B'}, [1; 2], 'VariableNames', {'Name', 'Value'});
uit = uitable(fig, 'Data', tdata, 'Position', [250 50 200 150]);

% tab group
tg = uitabgroup(fig);
t1 = uitab(tg, 'Title', 'Tab 1');
t2 = uitab(tg, 'Title', 'Tab 2');
plot(uiaxes(t1), 1:10);

% input dialog
answer = inputdlg({'Name:', 'Age:'}, 'User Info', [1 50; 1 10]);
if ~isempty(answer)
    name = answer{1};
    age = str2double(answer{2});
end

% question dialog
btn = questdlg('Save changes?', 'Confirm', 'Yes', 'No', 'Cancel', 'Yes');
switch btn
    case 'Yes', save_data();
    case 'No',  % discard
    case 'Cancel', return;
end

% list dialog
[sel, ok] = listdlg('PromptString', 'Select items:', ...
                    'ListString', {'Apple', 'Banana', 'Cherry'}, ...
                    'SelectionMode', 'multiple');

% file dialogs
[file, path] = uigetfile('*.csv', 'Select CSV');
[file, path] = uiputfile('*.mat', 'Save as');
fullpath = fullfile(path, file);

% message dialog
msgbox('Done!', 'Info', 'help');
errordlg('Something went wrong', 'Error');
warndlg('Check this', 'Warning');

% Run as standalone app
%   1. In App Designer: Designer Tab > App Store Details > Save as .mlapp
%   2. Share the .mlapp file — others open it in MATLAB

% Create standalone app (requires MATLAB Compiler)
%   >> mcr.build('myApp.mlapp')        % R2024a+
%   or
%   >> deploytool                       % legacy

% Web App (requires MATLAB Web App Server)
%   >> compiler.build.webApp('myApp.mlapp')

% Export to code (readable .m)
%   In App Designer: Designer Tab > Export > Export to .m file

% Pack as toolbox (for sharing with colleagues)
%   >> matlab.addons.toolbox.packageToolbox('myApp.prj)

% Tips:
%   - Use try/catch in callbacks to avoid silent failures
%   - Set app.UIFigure.Name to identify your app
%   - Test on a clean MATLAB session before deploying

fid = fopen('data.txt', 'r');   % r, w, a, r+, w+
if fid == -1; error('Cannot open file'); end

% read line by line
while ~feof(fid)
  line = fgetl(fid);              % without newline
  fprintf('%s\n', line);
end

% formatted read
% data.txt: "Alice 30 90.5"
A = fscanf(fid, '%s %d %f', [3 1]);

% formatted write
fprintf(fid, '%s, %d, %.2f\n', 'Bob', 25, 85.3);

% binary I/O
data = fread(fid, [10 5], 'double');  % 10x5 doubles
fwrite(fid, data, 'double');

% position
fseek(fid, 0, 'bof');             % beginning
fseek(fid, 0, 'eof');             % end
pos = ftell(fid);

fclose(fid);                      % ALWAYS close

fid = fopen('data.csv', 'r');

% skip header line
header = fgetl(fid);

% parse remaining
% Format: Alice,30,90.5,B+
C = textscan(fid, '%s %d %f %s', ...
             'Delimiter', ',', ...
             'MultipleDelimsAsOne', true, ...
             'TreatAsEmpty', {'NA', 'N/A'}, ...
             'HeaderLines', 0);
fclose(fid);

names = C{1}; ages = C{2}; scores = C{3}; grades = C{4};

% alternative: read whole file at once
str = fileread('data.txt');

% from string
C = textscan('1 2 3\n4 5 6', '%d %d %d');
% C{1} = [1; 4], C{2} = [2; 5], C{3} = [3; 6]

% collectoutput (group columns of same type)
C = textscan(fid, '%s %d %d %f', 'CollectOutput', true);
% C{2} = [ages, more_ints] as one matrix

A = magic(5);
B = struct('name', 'Alice', 'age', 30);
C = rand(100);

% save variables
save('data.mat', 'A', 'B', 'C');
save('data.mat', 'A', 'B', '-v7.3');  % large files > 2GB
save('data.mat', '-append', 'C');     % add to existing

% load
loaded = load('data.mat');       % struct with variables
loaded.A
loaded.B.age

% or load into workspace
load('data.mat');                % A, B, C now in workspace

% partial load (huge files)
info = whos('-file', 'bigdata.mat');
loaded = load('bigdata.mat', 'A');  % only A

% save specific format
save('data.txt', 'A', '-ascii', '-double');
save('data.csv', 'A', '-ascii', '-delimiter', ',');

% compressed
save('data.mat', 'A', '-v7.3', '-nocompression');

% HDF5 (great for large scientific datasets)
% write
h5create('data.h5', '/dataset1', [100 50 20]);
data = rand(100, 50, 20);
h5write('data.h5', '/dataset1', data);

% read
info = h5info('data.h5');
data = h5read('data.h5', '/dataset1');
partial = h5read('data.h5', '/dataset1', [1 1 1], [10 10 5]);

% append (extendable datasets)
h5create('data.h5', '/growable', [100 Inf], 'ChunkSize', [100 100]);
for k = 1:5
  h5write('data.h5', '/growable', rand(100, 100), [1 (k-1)*100+1]);
end

% NetCDF (climate/ocean data)
ncid = netcdf.open('data.nc', 'NC_NOWRITE');
varid = netcdf.inqVarID(ncid, 'temperature');
data = netcdf.getVar(ncid, varid);
netcdf.close(ncid);

% FITS (astronomy)
fitsdisp('image.fits');
data = fitsread('image.fits');

% JSON (R2016b+)
data = struct('name', 'Alice', 'age', 30, 'scores', [90 80 85]);
txt = jsonencode(data, 'PrettyPrint', true);
disp(txt);

decoded = jsondecode(txt);
decoded.name                     % 'Alice'

% JSON with cell arrays for mixed types
data = struct('items', {{'apple', 42, [1 2 3]}});
txt = jsonencode(data);

% read JSON file
txt = fileread('config.json');
config = jsondecode(txt);

% XML
doc = xmlread('data.xml');
root = doc.getDocumentElement;
children = root.getChildNodes;
for k = 0:children.getLength - 1
  node = children.item(k);
  if node.getNodeType == doc.ELEMENT_NODE
    fprintf('%s: %s\n', node.getNodeName, char(node.getTextContent));
  end
end

% write XML
doc = com.mathworks.xml.XMLUtils.createDocument('root');
root = doc.getDocumentElement;
child = doc.createElement('item');
child.appendChild(doc.createTextNode('hello'));
root.appendChild(child);
xmlwrite('out.xml', doc);

% requires Parallel Computing Toolbox
% start pool
pool = gcp('nocreate');
if isempty(pool), parpool; end

% parfor: parallel for loop
N = 1000;
results = zeros(1, N);
parfor k = 1:N
  results(k) = some_expensive_function(k);
end

% constraints:
%   - iterations must be independent (no order dependency)
%   - body cannot contain break/return
%   - variables classified as sliced, broadcast, reduction, temp

% reduction variables
total = 0;
parfor k = 1:N
  total = total + compute(k);
end

% sliced variables (independent indexing)
data = rand(N, 100);
out = zeros(N, 1);
parfor k = 1:N
  out(k) = mean(data(k, :));
end

% spmd: single program, multiple data
spmd
  % code runs on every worker
  % labindex: this worker's ID (1 to numlabs)
  % numlabs: total workers
  fprintf('Worker %d of %d\n', labindex, numlabs);

  % each worker computes a chunk
  local_data = rand(100, 1);
  local_sum = sum(local_data);
end

% combine results
total = sum([local_sum{:}]);

% distributed arrays (split across workers)
D = distributed.rand(10000, 10000);  % huge matrix
local_size = size(D, 'local')        % each worker's chunk
total_size = size(D)                  % logical total

% operations auto-distribute
S = D * D';
S_local = gather(S);                  % pull to client (if fits)

% Composite (cell-like, one per worker)
C = Composite();
spmd
  C{labindex} = magic(labindex + 1);
end
C{1}                              % access worker 1's value

% requires GPU support (most NVIDIA GPUs)
% move data to GPU
A = gpuArray(rand(10000));
B = gpuArray(rand(10000));

% operations run on GPU
C = A * B;                        % still on GPU
D = sin(A) + cos(B);
E = sum(A, 1);

% bring back to CPU
C_cpu = gather(C);

% check GPU
gpuDevice                        % info about current GPU
gpuDeviceCount                   % number of available GPUs
gpuDevice(1)                     % select device 1

% element-wise ops are typically fastest on GPU
% matrix multiply also very fast (cuBLAS)

% custom kernels (PTX)
k = parallel.gpu.CUDAKernel('my_kernel.ptx', 'my_kernel.cu');
result = feval(k, A, B);

% arrayfun on GPU (element-wise custom functions)
f = @(x) x^2 + sin(x);
result = arrayfun(f, A);

% run a function in background
job = batch(@my_function, 0, {arg1, arg2});

% check status
job.State                        % 'queued', 'running', 'finished'
wait(job);                       % block until done

% get results
results = fetchOutputs(job);
delete(job);                     % clean up

% batch script
job = batch('my_script.m');

% with pool
job = batch(@my_func, 1, {x}, 'Pool', 4);  % 4 extra workers

% batch with attached files
job = batch(@my_func, 1, {x}, ...
            'AttachedFiles', {'data.mat', 'helper.m'}, ...
            'CurrentFolder', '/path/to/work');

% list jobs
jobs = findJob(pool);
for j = jobs
  fprintf('Job %d: %s\n', j.ID, j.State);
end

% parallel pool settings
parpool('local', 4);              % 4 workers
delete(gcp('nocreate'));          % shut down pool

% 1. Preallocate arrays
% BAD
for k = 1:1000
  x(k) = k^2;                    % grows array each iteration
end

% GOOD
x = zeros(1, 1000);
for k = 1:1000
  x(k) = k^2;
end

% 2. Vectorize instead of looping
% BAD
for k = 1:length(x)
  y(k) = sin(x(k));
end

% GOOD
y = sin(x);                      % one vectorized call

% 3. Use JIT-friendly patterns
% simple loops are now fast (JIT), but vectorize for clarity

% 4. Profile to find bottlenecks
profile on;
my_function();
profile viewer;

% 5. Prefer single precision for huge data
A = rand(10000, 'single');       % half the memory

% 6. MEX for critical sections
% compile C code: mex my_func.c

% file: Point.m
classdef Point
  properties
    x = 0
    y = 0
  end

  properties (Dependent)
    r                              % computed on access
  end

  methods
    function obj = Point(x, y)
      if nargin > 0
        obj.x = x;
        obj.y = y;
      end
    end

    function r = get.r(obj)
      r = sqrt(obj.x^2 + obj.y^2);
    end

    function obj = move(obj, dx, dy)
      obj.x = obj.x + dx;
      obj.y = obj.y + dy;
    end

    function disp(obj)
      fprintf('Point(%.2f, %.2f), r=%.2f\n', obj.x, obj.y, obj.r);
    end
  end

  methods (Static)
    function p = origin()
      p = Point(0, 0);
    end
  end
end

% Value class (default) — copied on assignment
classdef VPoint
  properties; x; end
  methods
    function obj = VPoint(x); obj.x = x; end
    function obj = setX(obj, x); obj.x = x; end
  end
end

p = VPoint(1);
p2 = p;                           % copy
p2 = setX(p2, 5);
disp(p.x)                         % still 1 (unchanged)

% Handle class — passed by reference
classdef HPoint < handle
  properties; x; end
  methods
    function obj = HPoint(x); obj.x = x; end
    function setX(obj, x); obj.x = x; end   % no return needed
  end
end

h = HPoint(1);
h2 = h;                           % same object
h2.setX(5);
disp(h.x)                         % 5 (shared)

% handle class features:
%   - events and listeners
%   - destructor (delete method)
%   - copy() method

% base class
classdef Animal
  properties; name; end
  methods
    function obj = Animal(name); obj.name = name; end
    function speak(obj)
      error('Abstract: subclass must override');
    end
    function describe(obj)
      fprintf('%s says: ', obj.name);
      obj.speak();
      fprintf('\n');
    end
  end
end

% subclass
classdef Dog < Animal
  methods
    function obj = Dog(name); obj = obj@Animal(name); end
    function speak(obj); fprintf('Woof'); end
  end
end

classdef Cat < Animal
  methods
    function obj = Cat(name); obj = obj@Animal(name); end
    function speak(obj); fprintf('Meow'); end
  end
end

% polymorphism
animals = {Dog('Rex'), Cat('Whiskers'), Dog('Buddy')};
for k = 1:length(animals)
  animals{k}.describe();
end

% multiple inheritance
classdef FlyingDog < Dog & IFlyable
  ...
end

classdef Thermometer < handle
  properties
    temperature = 20
  end

  events
    temperatureChanged
    overheat
  end

  methods
    function set.temperature(obj, t)
      obj.temperature = t;
      notify(obj, 'temperatureChanged');
      if t > 100
        notify(obj, 'overheat');
      end
    end

    function obj = Thermometer()
      % add listener to own event
      addlistener(obj, 'overheat', @(s,e) disp('WARNING: too hot!'));
    end
  end
end

% usage
t = Thermometer();
lh = addlistener(t, 'temperatureChanged', @(s,e) ...
  fprintf('Now: %.1f\n', s.temperature));
lh2 = addlistener(t, 'overheat', @onOverheat);

t.temperature = 25;               % triggers temperatureChanged
t.temperature = 105;              % triggers both

% custom event data
classdef OverheatData < event.EventData
  properties; temp; end
end
% notify(obj, 'overheat', OverheatData(t));

% enumeration class
classdef Color
  enumeration
    Red [1 0 0]
    Green [0 1 0]
    Blue [0 0 1]
  end
  properties
    rgb
  end
  methods
    function obj = Color(rgb)
      obj.rgb = rgb;
    end
  end
end

c = Color.Red;
disp(c.rgb)                       % [1 0 0]

% switch on enum
switch c
  case Color.Red; disp('red');
  case Color.Green; disp('green');
end

% property attributes
classdef Account
  properties (SetAccess = private)
    balance = 0                   % readable outside, settable only inside
  end
  properties (Access = {?Account, ?Bank})
    internal_id                   % only Account and Bank classes
  end
  properties (Constant)
    PI = 3.14159                  % compile-time constant
  end
  properties (Hidden)
    cache                         % not shown by disp/struct
  end
end