EC23S01 NUMERICAL AND SIGNAL PROCESSING PRACTICE THROUGH
PYTHON LABORATORY
Getting Started with Python and using the Plot command
Continuous-Time Signal
(a) Plot a sine wave x(t)=sin(2πft)x(t) = \sin(2\pi f t)x(t)=sin(2πft)
import numpy as np
import [Link] as plt
from scipy import signal
# Create a time array from 0 to 1 second, with 1000 samples
t = [Link](0, 1, 1000, endpoint=False)
# Create a sine wave signal
freq = 5 # 5 Hz
x = [Link](2 * [Link] * freq * t)
# Plot the signal
[Link](t, x)
[Link]("Sine Wave -5 Hz")
[Link]("Time [s]")
[Link]("Amplitude")
[Link](True)
[Link]()
(b) Square Wave
import numpy as np
import [Link] as plt
x_square = [Link](2 * [Link] * freq * t)
[Link](t, x_square)
[Link]("Square Wave -5 Hz")
[Link]("Time [s]")
[Link]("Amplitude")
[Link](True)
[Link]()
(c) Unit Step Function
import numpy as np
import [Link] as plt
from scipy import signal
x_square = [Link](2 * [Link] * freq * t)
t = [Link](-1, 1, 1000)
x = [Link](t, 1)
[Link](t, x)
[Link]("Unit Step Function")
1. Scaling a Signal
To multiply each sample of a signal x[n]x[n]x[n] by a constant aaa:
y[n]=a⋅x[n]y[n
2. Signal Addition
Adding two signals x[n] and y[n]
z[n]=x[n]+y[n]
3. Reversing a Signal
To reverse a signal x[n]→x[−n]
xrev[n]=x[−n] (represented by reversing the list)
4. Shifting a Signal
Right Shift by k: y[n]=x[n−k]
Left Shift by k: y[n]=x[n+k]
Procedure:
Create a time vector n using range() or a list.
Define discrete-time signals (step, ramp) using list comprehensions.
Apply signal operations: scaling, addition, reversal, shifting.
Use [Link]() to plot signals.
Interpret the behavior of each signal and operation.
***************
import [Link] as plt
# 1. Creating basic signals using lists
n = list(range(0, 10)) # Time indices
# Unit step signal
unit_step = [1for i in n]
# Ramp signal
ramp = [i for i in n]
# 2. Scaling a signal
scaled_ramp = [2* x for x in ramp]
# 3. Adding two signals
combined_signal = [u + r for u, r in zip(unit_step, ramp)]
# 4. Reversing a signal
reversed_ramp = ramp[::-1]
# 5. Shifting a signal (Right shift by 2)
right_shifted = [0, 0] + ramp[:-2]
# Plotting all signals
[Link](figsize=(12, 8))
[Link](3, 2, 1)
[Link](n, unit_step)
[Link]("Unit Step Signal")
[Link]("n"); [Link]("Amplitude")
[Link](3, 2, 2)
[Link](n, ramp)
[Link]("Ramp Signal")
[Link](3, 2, 3)
[Link](n, scaled_ramp)
[Link]("Scaled Ramp Signal")
[Link](3, 2, 4)
[Link](n, combined_signal)
[Link]("Combined Signal (Step + Ramp)")
[Link](3, 2, 5)
[Link](n, reversed_ramp)
[Link]("Reversed Ramp Signal")
[Link](3, 2, 6)
[Link](n, right_shifted)
[Link]("Right Shifted Ramp Signal")
plt.tight_layout()
[Link]()
2. Load data from files and Plotting data
(a) sine wave import
(b) Square Wave
import numpy as np
import [Link] as plt
from scipy import signal
# Time vector
t = [Link](0, 1, 500, endpoint=False)
# Generate square wave
square_wave = [Link](2 * [Link] * 5 * t) # 5 Hz frequency
# Plotting
[Link](t, square_wave)
[Link]("5 Hz Square Wave")
[Link]("Time [s]")
[Link]("Amplitude")
[Link](True)
[Link]()
(c) Unit Step Function
import numpy as np
import [Link] as plt
# Time range
t = [Link](-5, 5, 1000)
# Unit step function
u = [Link](t >= 0, 1, 0)
# Plot
[Link](figsize=(8, 4))
[Link](t, u, label='Unit Step Function u(t)')
[Link]('Time')
[Link]('Amplitude')
[Link]('Unit Step Function')
[Link](True)
[Link]()
[Link]()
Exp. No.: 3 Getting Started with Date:
Lists
import math
import [Link] as plt
# Define n as a list
n = list(range(0, 10))
# Rectangular pulse x1[n] = 1 for n=0 to 4
x1 = list([1] * 5)
# Exponential decay signal x2[n] = exp(-0.5 * n)
x2 = [[Link](-0.5 * i) for i in n]
# Perform linear convolution manually using lists
m = len(x1)
n2 = len(x2)
conv_result = list([0] * (m + n2 - 1))
for i in range(m + n2 - 1):
for j in range(m):
if 0 <= i - j < n2:
conv_result[i] += x1[j] * x2[i - j]
# Plotting
[Link](figsize=(10, 6))
# x1[n]
[Link](3, 1, 1)
[Link](range(len(x1)), x1, linefmt='b-', markerfmt='bo', basefmt=" ", label='Rectangular Pulse')
[Link]()
# x2[n]
[Link](3, 1, 2)
[Link](range(len(x2)), x2, linefmt='g-', markerfmt='go', basefmt=" ", label='Exponential Decay')
[Link]()
# Convolution result
[Link](3, 1, 3)
[Link](range(len(conv_result)), conv_result, linefmt='r-', markerfmt='ro', basefmt="",
label='Convolution Result')
[Link]()
[Link]('n')
[Link](True)
plt.tight_layout()
[Link]()
***linear conv
import [Link] as plt
# Generate Input Signals using list()
x = list([1, 2, 3, 4])
h = list([1, 1, 1, -1])
# Lengths
m = len(x)
n = len(h)
# Zero Padding using list() and +
X = x + list([0]*n)
H = h + list([0]*m)
# Convolution Operation using lists
Y = list([0]*(m + n - 1))
for i in range(m + n - 1):
for j in range(m):
if 0 <= i - j < n + m:
Y[i] += X[j] * H[i - j]
# Plot results
[Link](figsize=(10, 6))
[Link](3, 1, 1)
[Link](range(len(x)), x, linefmt='b-', markerfmt='b^',
basefmt='k')
[Link]('n')
[Link]('x[n]')
[Link](True)
[Link](3, 1, 2)
[Link](range(len(h)), h, linefmt='m-', markerfmt='ms',
basefmt='k')
[Link]('n')
[Link]('h[n]')
[Link](True)
[Link](3, 1, 3)
[Link](range(len(Y)), Y, linefmt='r-', markerfmt='ro',
basefmt='k')
[Link]('n')
[Link]('Y[n]')
[Link](True)
plt.tight_layout()
[Link]()
circular
************
import numpy as np
import [Link] as plt
# Input sequences using list()
x = list([1, 2, 3, 4])
h = list([1, 1, 1, -1])
# Length for circular convolution = max(len(x), len(h))
N = max(len(x), len(h))
# Zero pad both sequences to length N using list() and +
x_padded = x + list([0] * (N - len(x)))
h_padded = h + list([0] * (N - len(h)))
# Perform FFT-based circular convolution
X = [Link](x_padded)
H = [Link](h_padded)
Y = [Link](X * H).real # Keep real part
# Convert result to list
Y = list(Y)
# Plot results
[Link](figsize=(10, 6))
[Link](3, 1, 1)
[Link](range(len(x)), x, linefmt='b-', markerfmt='b^',
basefmt='k')
[Link]('n')
[Link]('x[n]')
[Link](True)
[Link](3, 1, 2)
[Link](range(len(h)), h, linefmt='m-', markerfmt='ms',
basefmt='k')
[Link]('n')
[Link]('h[n]')
[Link](True)
[Link](3, 1, 3)
[Link](range(len(Y)), Y, linefmt='r-', markerfmt='ro',
basefmt='k')
[Link]('n')
[Link]('Y[n]')
[Link](True)
plt.tight_layout()
[Link]()
import numpy as np
import [Link] as plt
from [Link] import butter, filtfilt
# Sample input signal
x = [1, 2, 3, 4] # You can use any real or complex signal
N = len(x)
X = [Link](x) # Use numpy's FFT for efficiency and correctness
# Print the DFT result
print("DFT of x[n]:")
for i, val in enumerate(X):
print(f"X[{i}] = {val}")
# Parameters
fs_list = [10, 20, 50]
f, duration, noise_std = 5, 2, 0.5
cutoff = 3
# Define low-pass filter
def lowpass_filter(data, cutoff, fs):
b, a = butter(4, cutoff / (0.5 * fs), btype='low')
return filtfilt(b, a, data)
# High-resolution time axis
time_fine = [Link](0, duration, 1000)
clean_signal = [Link](2 * [Link] * f * time_fine)
# Plotting
[Link](figsize=(10, 8))
for i, fs in enumerate(fs_list):
# Coarse sampling for noisy signal
time = [Link](0, duration, int(fs * duration), endpoint=False)
noisy_signal = [Link](2 * [Link] * f * time) + [Link](0, noise_std, len(time))
filtered_signal = lowpass_filter(noisy_signal, cutoff, fs)
# Interpolate to high-resolution time base
noisy_interp = [Link](time_fine, time, noisy_signal)
filtered_interp = [Link](time_fine, time, filtered_signal)
# Plot
[Link](len(fs_list), 1, i + 1)
[Link](time_fine, noisy_interp, 'r--', alpha=0.5, label='Noisy')
[Link](time_fine, filtered_interp, 'b', label='Filtered')
[Link](time_fine, clean_signal, 'k', alpha=0.6, label='Clean')
[Link](f"Sampling Rate: {fs} Hz")
[Link]("Time (s)")
[Link]("Amplitude")
[Link](loc="upper right")
plt.tight_layout()
[Link]()
Exp. No.: 4 Getting started with for, Date:
If, While loops
Fourier Transform of Basic Signals
Write a Python program to compute the Discrete Fourier Transform (DFT) of a sinusoidal signal
x(t)=sin(2π10t). Use [Link] to perform the Fourier Transform and plot both the time-domain
signal and the magnitude spectrum. Discuss how the frequency content is reflected in the spectrum.
DFT of the signal using inbuilt function
import numpy as np
import [Link] as plt
t = [Link](0, 1,500)
x = [Link](2 * [Link] * 100 * t)+[Link](2 * [Link] * 200 * t)
X_f = [Link](x)
frequencies = [Link](len(t), d=t[1] - t[0])
#print(frequencies)
magnitude_spectrum = [Link](X_f)
#print(magnitude_spectrum)
[Link](figsize=(12, 6))
[Link](2, 1, 1)
[Link](t, x, label='x(t) = sin(2π10t)')
[Link]('Time (s)')
[Link]('Amplitude')
[Link]('Time-Domain Signal')
[Link](True)
[Link](2, 1, 2)
[Link](frequencies[:len(frequencies)//2],
magnitude_spectrum[:len(frequencies)//2])
[Link]('Frequency (Hz)')
[Link]('Magnitude')
[Link]('Magnitude Spectrum')
[Link](True)
plt.tight_layout()
[Link]()
Dft of the signal using formalue
x = [1, 2, 3, 4] # You can use any real or complex signal
N = len(x)
X = []
# DFT formula:
#X[k] = sum_n=0_to_N-1 x[n] * exp(-j*2*pi*k*n/N)
for k in range(N):
real = 0
imag = 0
for n in range(N):
angle = 2 * [Link] * k * n / N
real += x[n] * [Link](-angle)
imag += x[n] * [Link](-angle)
[Link](complex(real, imag))
#Print the DFT result
print("DFT of x[n]:")
for i, val in enumerate(X):
print(f"X[{i}] = {val}")
X[0] = (10+0j)
X[1] = (-2.0000000000000004+1.9999999999999996j)
X[2] = (-2-9.797174393178826e-16j)
X[3] = (-1.9999999999999982-2.000000000000001j)
Using for statement signal filtering
**
import numpy as np
import [Link] as plt
from [Link] import butter, filtfilt
# Parameters
fs_list = [10, 20, 50]
f, duration, noise_std = 5, 2, 0.5
cutoff = 3
# Define low-pass filter
def lowpass_filter(data, cutoff, fs):
b, a = butter(4, cutoff / (0.5 * fs), btype='low')
return filtfilt(b, a, data)
# High-resolution time axis
time_fine = [Link](0, duration, 1000)
clean_signal = [Link](2 * [Link] * f * time_fine)
# Plotting
[Link](figsize=(10, 8))
for i, fs in enumerate(fs_list):
# Coarse sampling for noisy signal
time = [Link](0, duration, int(fs * duration), endpoint=False)
noisy_signal = [Link](2 * [Link] * f * time) + [Link](0,
noise_std, len(time))
filtered_signal = lowpass_filter(noisy_signal, cutoff, fs)
# Interpolate to high-resolution time base
noisy_interp = [Link](time_fine, time, noisy_signal)
filtered_interp = [Link](time_fine, time, filtered_signal)
# Plot
[Link](len (fs_list), 1, i+1)
print(i)
print(fs_list)
[Link](time_fine, noisy_interp, 'r--', alpha=0.5, label='Noisy')
[Link](time_fine, filtered_interp, 'b', label='Filtered')
[Link](time_fine, clean_signal, 'k', alpha=0.6, label='Clean')
[Link](f"Sampling Rate: {fs} Hz")
[Link]("Time (s)")
[Link]("Amplitude")
[Link](loc="upper right")
plt.tight_layout()
[Link]()
# Original signal (can be any list of values)
signal = [1, 3, 7, 12, 5, -4, -10, 2]
# Clipping thresholds
min_val = -5
max_val = 10
# Clipped signal (using 'if')
clipped_signal = []
for sample in signal:
if sample > max_val:
clipped_signal.append(max_val)
elif sample < min_val:
clipped_signal.append(min_val)
else:
clipped_signal.append(sample)
print("Original Signal: ", signal)
print("Clipped Signal: ", clipped_signal)
while loop for signal normalizxation
signal = [4, 8, 15, 16, 23, 42]
# Find min and max
min_val = min(signal)
max_val = max(signal)
# Normalized signal (initialize)
normalized_signal = []
# Counter for while loop
i = 0
while i<len(signal):
# Apply normalization formula
norm = (signal[i] - min_val) / (max_val - min_val)
normalized_signal.append(norm)
i += 1
# Display results
print("Original Signal: ", signal)
print("Normalized Signal: ", normalized_signal)
Original Signal: [4, 8, 15, 16, 23, 42]
Normalized Signal: [0.0, 0.10526315789473684, 0.2894736842105263,
0.3157894736842105, 0.5, 1.0]
&&&&&&&&&&&&&&&&&&&&&&&&&
5/Getting started with files and arrays Aliasing in Real-World Signals
import numpy as np
import [Link] as plt
from [Link] import resample
# Original signal parameters
sr = 44100
t = [Link](0, 1, sr, endpoint=False)
original_signal = [Link](2 * [Link] * 440 * t) + 0.5 * [Link](2 * [Link] *
880 * t)
# Downsampling rates
downsampling_rates = [8000, 16000, 44100]
# Create figure
[Link](figsize=(12, 8))
for idx, rate in enumerate(downsampling_rates):
# Downsample the signal
downsampled_signal = resample(original_signal, rate)
# Compute FFT and frequency axis
frequencies = [Link](len(downsampled_signal), d=1/rate)
Y_f = [Link](downsampled_signal)
# Plot time-domain signal
[Link](3, 2, 2 * idx + 1)
[Link](downsampled_signal)
[Link](f'Downsampled Signal at {rate} Hz')
[Link]('Sample')
[Link]('Amplitude')
# Plot frequency-domain spectrum
[Link](3, 2, 2 * idx + 2)
[Link](frequencies[:len(frequencies)//2], [Link](Y_f)[:len(Y_f)//2])
[Link](f'Frequency Spectrum at {rate} Hz')
[Link]('Frequency (Hz)')
[Link]('Magnitude')
# Final layout adjustment
plt.tight_layout()
[Link]()
CODE: Using librosa
import librosa
import numpy as np
import [Link] as plt
from [Link] import Audio
# Load audio file (example trumpet sound)
signal, sr = [Link]([Link]('trumpet'))
# Function to plot frequency spectrum
def plot_spec(sig, sr, tle="Spectrum"):
D = [Link]([Link](sig))
freqs = [Link](len(D), 1/sr)
[Link](figsize=(8, 4))
[Link](freqs[:len(freqs)//2], D[:len(D)//2])
[Link](tle) # Fixed typo from 'tle'
[Link]("Frequency (Hz)")
[Link]("Amplitude")
[Link](True)
[Link]()
# Plot spectrum of original signal
plot_spec(signal, sr, "Original Signal Spectrum")
# Downsample rates
rates = [8000, 16000, 44100]
# Downsample, plot spectrum, and play audio
for r in rates:
ds_signal = [Link](y=signal, orig_sr=sr, target_sr=r)
plot_spec(ds_signal, r, f"Downsampled to {r} Hz")
import librosa
import numpy as np
import [Link] as plt
from [Link] import Audio
# Load audio file (example trumpet sound)
signal, sr = [Link]([Link]('trumpet'))
# Function to plot frequency spectrum
def plot_spec(sig, sr, tle="Spectrum"):
D = [Link]([Link](sig))
freqs = [Link](len(D), 1/sr)
[Link](figsize=(8, 4))
[Link](freqs[:len(freqs)//2], D[:len(D)//2])
[Link](tle) # Fixed typo from 'tle'
[Link]("Frequency (Hz)")
[Link]("Amplitude")
[Link](True)
[Link]()
# Plot spectrum of original signal
plot_spec(signal, sr, "Original Signal Spectrum")
# Downsample rates
rates = [8000, 16000, 44100]
# Downsample, plot spectrum, and play audio
for r in rates:
ds_signal = [Link](y=signal, orig_sr=sr, target_sr=r)
plot_spec(ds_signal, r, f"Downsampled to {r} Hz")
Power Spectral Density (PSD) - Using Numpy Array
Write a Python program to compute the Power Spectral Density (PSD) of a noisy signal consisting of
a 5 Hz sinusoid and Gaussian noise. Use the matplotlib library to plot the PSD and explain how the
noise and signal components appear in the frequency
import numpy as np
import [Link] as plt
from [Link] import welch
# Sampling settings
fs = 100 # Sampling frequency (Hz)
t_duration = 2 # Duration in seconds
f_signal = 5 # Frequency of the sine wave (Hz)
# Time axis - creates a 1D NumPy array of time points
t = [Link](0, t_duration, 1/fs)
# Signal + noise
# Creates a sine wave using the time array t
signal = [Link](2 * [Link] * f_signal * t)
# Generates Gaussian noise as a NumPy array of the same length as t
noise = [Link](0, 0.5, len(t))
# Element-wise addition of two arrays (signal + noise).
noisy_signal = signal + noise
# Power Spectral Density using Welch method
frequencies, psd = welch(noisy_signal, fs, nperseg=256)
# Plot
[Link](figsize=(10, 6))
# Plotting two arrays: frequency values vs. power values.
[Link](frequencies, psd, label='PSD of Noisy Signal')
[Link]('Power Spectral Density (PSD)')
[Link]('Frequency (Hz)')
[Link]('Power/Frequency (dB/Hz)')
[Link](True)
[Link]()
[Link](0, 30)
[Link](1e-5, 1e1)
plt.tight_layout()
[Link]()
import numpy as np
import [Link] as plt
from [Link] import welch
# Sampling settings
fs = 100 # Sampling frequency (Hz)
t_duration = 2 # Duration in seconds
f_signal = 5 # Frequency of the sine wave (Hz)
# Time axis - creates a 1D NumPy array of time points
t = [Link](0, t_duration, 1/fs)
# Signal + noise
# Creates a sine wave using the time array t
signal = [Link](2 * [Link] * f_signal * t)
# Generates Gaussian noise as a NumPy array of the same length as t
noise = [Link](0, 0.5, len(t))
# Element-wise addition of two arrays (signal + noise).
noisy_signal = signal + noise
# Power Spectral Density using Welch method
frequencies, psd = welch(noisy_signal, fs, nperseg=256)
# Plot
[Link](figsize=(10, 6))
# Plotting two arrays: frequency values vs. power values.
[Link](frequencies, psd, label='PSD of Noisy Signal')
[Link]('Power Spectral Density (PSD)')
[Link]('Frequency (Hz)')
[Link]('Power/Frequency (dB/Hz)')
[Link](True)
[Link]()
[Link](0, 30)
[Link](1e-5, 1e1)
plt.tight_layout()
[Link]()
Statistics using Python
import numpy as np
import [Link] as plt
from [Link] import skew, kurtosis
# Generate a noisy sine signal
fs = 1000
t = [Link](0, 1, fs, endpoint=False)
signal = [Link](2 * [Link] * 5 * t)
noise = [Link](0, 0.5, size=fs)
noisy_signal = signal + noise
# --- Statistics ---
mean_val = [Link](noisy_signal)
std_val = [Link](noisy_signal)
rms_val = [Link]([Link](noisy_signal**2))
skew_val = skew(noisy_signal)
kurt_val = kurtosis(noisy_signal)
# --- Output Statistics ---
print(f"Mean: {mean_val:.4f}")
print(f"Standard Deviation: {std_val:.4f}")
print(f"RMS: {rms_val:.4f}")
print(f"Skewness: {skew_val:.4f}")
print(f"Kurtosis: {kurt_val:.4f}")
# --- Plot Signal & Histogram ---
[Link](figsize=(12, 5))
[Link](1, 2, 1)
[Link](t, noisy_signal)
[Link]("Noisy Signal")
[Link]("Time (s)")
[Link]("Amplitude")
[Link](1, 2, 2)
[Link](noisy_signal, bins=50, density=True, color='skyblue', edgecolor='black')
[Link]("Signal Amplitude Histogram")
[Link]("Amplitude")
[Link]("Probability Density")
plt.tight_layout()
[Link]()
Linear combination of vectors and Computation of determinant, rank of a matrix
import numpy as np
# Define vectors
v1 = [Link]([1, 2])
v2 = [Link]([3, 4])
c1 = 2
c2 = -1
# Linear combination: v = c1*v1 + c2*v2
v = c1 * v1 + c2 * v2
print("Vector v1:", v1)
print("Vector v2:", v2)
print("Linear combination (2*v1 - 1*v2):", v)
#2[1 2 ] -1[3 4]
#2-3 4-4
#[-1 0]
import numpy as np
# Define matrix
A = [Link]([[2, 3, 1],
[4, 7, 7],
[6, 18, 22]])
#2(154-126)-3(88-42)+1(72-42)
#2(28) -138+30
#56-138+30
#-52
[1 2 3
142
265]
# Compute determinant
det = [Link](A)
# Compute rank
rank = [Link].matrix_rank(A)
print("Matrix A:\n", A)
print("Determinant of A:", round(det, 4))
print("Rank of A:", rank)
import numpy as np
def lu_decomposition(A):
n = len(A)
L = [Link]((n, n))
U = [Link]((n, n))
for i in range(n):
# Upper Triangular
for k in range(i, n):
sum_ = sum(L[i][j] * U[j][k] for j in range(i))
U[i][k] = A[i][k] - sum_
# Lower Triangular
for k in range(i, n):
if i == k:
L[i][i] = 1
else:
sum_ = sum(L[k][j] * U[j][i] for j in range(i))
L[k][i] = (A[k][i] - sum_) / U[i][i]
return L, U
def forward_substitution(L, b):
n = len(b)
y = [Link](n)
for i in range(n):
y[i] = b[i] - sum(L[i][j] * y[j] for j in range(i))
return y
def backward_substitution(U, y):
n = len(y)
x = [Link](n)
for i in reversed(range(n)):
x[i] = (y[i] - sum(U[i][j] * x[j] for j in range(i+1, n))) / U[i][i]
return x
# Main function
A = [Link]([[2, 3, 1],
[4, 7, 7],
[6, 18, 22]], dtype=float)
b = [Link]([1, 2, 3], dtype=float)
L, U = lu_decomposition(A)
y = forward_substitution(L, b)
x = backward_substitution(U, y)
print("Matrix L:\n", L)
print("Matrix U:\n", U)
print("Solution x:\n", x)
import numpy as np
def lu_decomposition(A):
n = len(A)
L = [Link]((n, n))
U = [Link]((n, n))
for i in range(n):
# Upper Triangular
for k in range(i, n):
sum_ = sum(L[i][j] * U[j][k] for j in range(i))
U[i][k] = A[i][k] - sum_
# Lower Triangular
for k in range(i, n):
if i == k:
L[i][i] = 1
else:
sum_ = sum(L[k][j] * U[j][i] for j in range(i))
L[k][i] = (A[k][i] - sum_) / U[i][i]
return L, U
def forward_substitution(L, b):
n = len(b)
y = [Link](n)
for i in range(n):
y[i] = b[i] - sum(L[i][j] * y[j] for j in range(i))
return y
def backward_substitution(U, y):
n = len(y)
x = [Link](n)
for i in reversed(range(n)):
x[i] = (y[i] - sum(U[i][j] * x[j] for j in range(i+1, n))) / U[i][i]
return x
# Main function
A = [Link]([[2, 3, 1], [4, 7, 7], [6, 18, 22]], dtype=float)
b = [Link]([1, 2, 3], dtype=float)
L, U = lu_decomposition(A)
y = forward_substitution(L, b)
x = backward_substitution(U, y)
print("Matrix L:\n", L)
print("Matrix U:\n", U)
print("Solution x:\n", x)
Generation of basic sequences using Python
import numpy as np
import [Link] as plt
n = [Link](-10, 11, 1)
# Unit Impulse
impulse = [Link](n == 0, 1, 0)
# Unit Step
step = [Link](n >= 0, 1, 0)
# Ramp
ramp = [Link](n >= 0, n, 0)
# Exponential
a = 0.8
expo = a ** n
# Sinusoidal
sin_seq = [Link](0.2 * [Link] * n)
# Plotting
[Link](figsize=(10, 10))
[Link](3, 2, 1)
[Link](n, impulse)
[Link]('Unit Impulse δ[n]')
[Link](3, 2, 2)
[Link](n, step)
[Link]('Unit Step u[n]')
[Link](3, 2, 3)
[Link](n, ramp)
[Link]('Ramp Sequence')
[Link](3, 2, 4)
[Link](n, expo)
[Link]('Exponential Sequence')
[Link](3, 2, 5)
[Link](n, sin_seq)
[Link]('Sinusoidal Sequence')
plt.tight_layout()
[Link]()
Generation of basic sequences using Python
import numpy as np
import [Link] as plt
n = [Link](-10, 11, 1)
# Unit Impulse
impulse = [Link](n == 0, 1, 0)
# Unit Step
step = [Link](n >= 0, 1, 0)
# Ramp
ramp = [Link](n >= 0, n, 0)
# Exponential
a = 0.8
expo = a ** n
# Sinusoidal
sin_seq = [Link](0.2 * [Link] * n)
# Plotting
[Link](figsize=(10, 10))
[Link](3, 2, 1)
[Link](n, impulse)
[Link]('Unit Impulse δ[n]')
[Link](3, 2, 2)
[Link](n, step)
[Link]('Unit Step u[n]')
[Link](3, 2, 3)
[Link](n, ramp)
[Link]('Ramp Sequence')
[Link](3, 2, 4)
[Link](n, expo)
[Link]('Exponential Sequence')
[Link](3, 2, 5)
[Link](n, sin_seq)
[Link]('Sinusoidal Sequence')
plt.tight_layout()
[Link]()
Spectral analysis of signals
import numpy as np
import [Link] as plt
# Generating a sample musical note signal
fs = 1100 # Sampling frequency (Hz)
duration = 2 # seconds
frequency = 440 # A4 note frequency (Hz)
t = [Link](0, duration, int(fs * duration), endpoint=False)
signal = [Link](2 * [Link] * frequency * t) + [Link](0, 1, len(t)) # Signal with noise
# Applying FFT
fft_result = [Link](signal)
freq = [Link]([Link][-1], d=1/fs)
# Plotting the spectrum
[Link](freq, [Link](fft_result))
[Link]('FFT of a Musical Note')
[Link]('Frequency (Hz)')
[Link]('Amplitude')
[Link]()
import numpy as np
import [Link] as plt
# Generating a sample musical note signal
fs = 1100 # Sampling frequency (Hz)
duration = 2 # seconds
frequency = 440 # A4 note frequency (Hz)
t = [Link](0, duration, int(fs * duration),
endpoint=False)
signal = [Link](2 * [Link] * frequency * t) +
[Link](0, 1, len(t)) # Signal with
noise
# Applying FFT
fft_result = [Link](signal)
freq = [Link]([Link][-1], d=1/fs)
# Plotting the spectrum
[Link](freq, [Link](fft_result))
[Link]('FFT of a Musical Note')
[Link]('Frequency (Hz)')
[Link]('Amplitude')
[Link]()