clear all; close all; clc;

%% Task-1

x = linspace(0.1, 5, 1000); % Define the range of x values
y = log(x); % Compute the function values
y0 = log(1); % Compute the value of the function at x=1

% Compute the Taylor series approximations
y2 = y0 + (x-1) - ((x-1).^2)/2;
y4 = y2 + (((x-1).^3)/3) - (((x-1).^4)/4);
y6 = y4 + (((x-1).^5)/5) - (((x-1).^6)/6);

% Plot the functions
figure;
plot(x, y, 'LineWidth', 2); % Plot the original function
hold on; % Hold the plot
plot(x, y2, '--', 'LineWidth', 2); % Plot the Taylor approximation with n=2
plot(x, y4, '-.', 'LineWidth', 2); % Plot the Taylor approximation with n=4
plot(x, y6, ':', 'LineWidth', 2); % Plot the Taylor approximation with n=6

% Add labels and a legend
xlabel('x');
ylabel('f(x)');
title('Taylor approximations of ln(x) around x=1');
legend('f(x)', 'n=2', 'n=4', 'n=6');


%% Task-2

x = linspace(0, 2*pi, 1000);
exact = sin(2*x);
approx = sin_approx(x);

figure;
plot(x, exact, 'LineWidth', 2);
hold on;
plot(x, approx, '--', 'LineWidth', 1);
xlabel('x');
ylabel('sin(2x)');
title('Exact and Approximate Values of sin(2x)');
legend('Exact', 'Approximation');

%% Task-3

n = 15000; % number of terms in the partial sum
s = zeros(n,1); % initialize the partial sum vector
x = 1:n;

for k = 2:n
    s(k) = s(k-1) + 1/(k*(log(k)^2));
end

figure;
plot(x, s);
title('Partial Sum of the Series \Sigma_{k=2}^{\infty} 1/(k(ln k)^2)')
xlabel('Number of Terms')
ylabel('Partial Sum')

%% Task-4

% Define the value of N
N = 11;

% Calculate the amplitude spectrum
n = 1:N;
An = abs((-1).^n + 1) ./ (pi*n);

% Calculate the phase spectrum
phi_n = zeros(1, N);
phi_n(1:2:end) = pi/2;

figure;

% Plot the amplitude spectrum
subplot(2,1,1)
stem(n, An, 'r', 'filled')
xlabel('n')
ylabel('Amplitude')
title('Amplitude Spectrum')

% Plot the phase spectrum
subplot(2,1,2)
stem(n, phi_n, 'b', 'filled')
xlabel('n')
ylabel('Phase (rad)')
title('Phase Spectrum')

% Reconstructed Waveform

T = 2;
N = 11;
x = linspace(-T/2, T/2, 1000);
f = 1/2;
for n = 1:N
    bn = (-((-1)^n)+1)/(pi*n);
    fn = bn*sin(n*pi*x);
    f = f + fn;
end

figure;
plot(x,f);
xlabel('x');
ylabel('f(x)');
title('Reconstructed Waveform');

%% Task-5

xn = [0,1,1,1,1,1,0,0]; % define the signal x[n]
X = fft(xn); % compute the DFT of x[n]

X_mag = abs(X); % compute the magnitude of X[k]
X_phase = angle(X); % compute the phase of X[k]

figure;
% plot the magnitude of X[k]
subplot(2,1,1);
stem(X_mag);
xlabel('k');
ylabel('|X[k]|');
title('Magnitude of DFT');

% plot the phase of X[k]
subplot(2,1,2);
stem(X_phase);
xlabel('k');
ylabel('Phase of X[k] (radians)');
title('Phase of DFT');

xi = ifft(X); % compute the IDFT of X[k]
n = 0:length(xi)-1; % define the sample indices

figure;
stem(n,xi);
xlabel('n');
ylabel('x[n]');
title('Inverse DFT');

%% Task-6

[x,Fs]= audioread(['GUI.m4a']);
L=size(x,1);

figure;
t=linspace(0,L/Fs,L);
subplot(2,1,1)
plot(t,x)
title("GUI audio signal w.r.t time")
xlabel("time (s)")
ylabel("amplitude (a.u.)")
Y = fft(x);
P2 = abs(Y/L);
P1 = P2(1:int32(L/2)+1);
P1(2:end-1) = 2*P1(2:end-1);
f = Fs*(0:(L/2))/L;
subplot(2,1,2)
plot(f,P1(2:end))
title("Single-Sided Amplitude Spectrum of GUI")
xlabel("f (Hz)")
ylabel("|P1(f)|")

%% Task-7

b = [0 1]; % numerator coefficients
a = [1 -1 -1]; % denominator coefficients

figure;
% Plot the poles and zeros
zplane(b,a);
hold on;

% Determine the ROC
z = roots(b);% finding zeros from the roots of b
p = roots(a); % finding poles from the roots of a
r = abs(p(end));
theta = 0:0.01:2*pi;
x = r*cos(theta);
y = r*sin(theta);

% Plot the ROC
patch(x,y,'red','FaceAlpha',0.1);

%% Task-8

A = [4 1 -7]; % A vector
B = [3 0.5 -7]; % B vector

% Dot product of A and B
dotAB = dot(A, B)

% Dot product of A and A
dotAA = dot(A, A)

% Dot product of B and A
dotBA = dot(B, A)

% Dot product of B and B
dotBB = dot(B, B)

% Cross product of A and B
crossAB = cross(A, B)

% Cross product of A and A
crossAA = cross(A, A)

% Cross product of B and A
crossBA = cross(B, A)

% Cross product of B and B
crossBB = cross(B, B)

%% Task-8 ii)

% (a)

A = [2 1 -3];
B = [1 0 -4];
C = [4 3 -1];

Ma = [A; B ; C];

rrefa = rref(Ma)

crossABdotC = dot(cross(A, B), C)


% (b)

A = [1 -3 2];
B = [2 -4 -1];
C = [3 2 -1];

Mb = [A; B ; C];

rrefb = rref(Mb)

crossABdotC2 = dot(cross(A, B), C)


%% Utility Function

function Y = sin_approx(X)
Y = zeros(size(X));

% sin_approx(x) returns an approximation of sin(2x) using Taylor series
% with a maximum error of 0.01.

for ii=1:numel(X)

    x = X(ii);

    % Set the maximum number of terms in the Taylor series
    n = 10;
    
    % Initialize the approximation to zero
    y = 0;
    
    % Compute the Taylor series
    for k = 0:n
        y = y + (-1)^k * (2*x)^(2*k+1) / factorial(2*k+1);
    end
    
    % Deal with the periodicity using mod function
    y = mod(y + pi, 2*pi) - pi;
    
    % Check if the maximum error is within the tolerance
    while abs(sin(2*x) - y) > 0.01
        n = n + 1;
        y = y + (-1)^n * (2*x)^(2*n+1) / factorial(2*n+1);
    end

    Y(ii) = y;

end

end