Signal Processing and Waveform Functions

def squarew(t, period=1):
    """
    Square wave function.
    
    Args:
        t: Time variable
        period: Period of the waveform (default: 1)
        
    Returns:
        Square wave value: +1 for first half of period, -1 for second half
    """

def trianglew(t, period=1):
    """
    Triangle wave function.
    
    Args:
        t: Time variable  
        period: Period of the waveform (default: 1)
        
    Returns:
        Triangle wave value: linear ramp from -1 to +1 and back
    """

def sawtoothw(t, period=1, rising=True):
    """
    Sawtooth wave function.
    
    Args:
        t: Time variable
        period: Period of the waveform (default: 1) 
        rising: True for rising sawtooth, False for falling (default: True)
        
    Returns:
        Sawtooth wave value: linear ramp within each period
    """

def unit_triangle(t, center=0, width=1):
    """
    Unit triangle function (triangular pulse).
    
    Args:
        t: Time variable
        center: Center position of triangle (default: 0)
        width: Base width of triangle (default: 1)
        
    Returns:
        Unit triangle: 1 at center, 0 outside width, linear in between
    """

def sigmoid(x):
    """
    Sigmoid (logistic) function.
    
    Args:
        x: Input value
        
    Returns:
        1 / (1 + exp(-x))
        Smooth S-curve from 0 to 1
    """

import mpmath
from mpmath import mp
import numpy as np  # For demonstration arrays

# Set precision
mp.dps = 25

# Example 1: Generate square wave
print("=== Square Wave ===")
t_values = [i * 0.1 for i in range(-20, 21)]  # -2 to 2 with step 0.1
square_values = [mp.squarew(t, period=1) for t in t_values]

print("Time\tSquare Wave")
for i in range(0, len(t_values), 4):  # Print every 4th value
    print(f"{t_values[i]:.1f}\t{square_values[i]}")

# Example 2: Triangle wave
print("\n=== Triangle Wave ===")
triangle_values = [mp.trianglew(t, period=2) for t in t_values]

print("Time\tTriangle Wave")
for i in range(0, len(t_values), 4):
    print(f"{t_values[i]:.1f}\t{triangle_values[i]}")

# Example 3: Sawtooth waves
print("\n=== Sawtooth Waves ===")
rising_sawtooth = [mp.sawtoothw(t, period=1.5, rising=True) for t in t_values]
falling_sawtooth = [mp.sawtoothw(t, period=1.5, rising=False) for t in t_values]

print("Time\tRising\tFalling")
for i in range(0, len(t_values), 4):
    print(f"{t_values[i]:.1f}\t{rising_sawtooth[i]:.3f}\t{falling_sawtooth[i]:.3f}")

# Example 4: Unit triangle function
print("\n=== Unit Triangle ===")
unit_tri_values = [mp.unit_triangle(t, center=0, width=2) for t in t_values]

print("Time\tUnit Triangle")
for i in range(len(t_values)):
    if unit_tri_values[i] > 0:  # Only print non-zero values
        print(f"{t_values[i]:.1f}\t{unit_tri_values[i]:.3f}")

# Example 5: Sigmoid function
print("\n=== Sigmoid Function ===")
x_values = [-3, -2, -1, 0, 1, 2, 3]
sigmoid_values = [mp.sigmoid(x) for x in x_values]

print("x\tSigmoid(x)")
for x, sig in zip(x_values, sigmoid_values):
    print(f"{x}\t{sig:.6f}")

# Example 6: Fourier series approximation using square wave
print("\n=== Fourier Series Approximation ===")
def fourier_square_wave(t, n_terms=5):
    """Approximate square wave using Fourier series."""
    result = 0
    for k in range(1, n_terms + 1):
        if k % 2 == 1:  # Only odd harmonics
            coefficient = 4 / (mp.pi * k)
            result += coefficient * mp.sin(2 * mp.pi * k * t)
    return result

t_test = [0, 0.25, 0.5, 0.75, 1.0]
print("Time\tSquare Wave\tFourier Approx (5 terms)")
for t in t_test:
    exact = mp.squarew(t)
    approx = fourier_square_wave(t, 5)
    print(f"{t:.2f}\t{exact:8.3f}\t{approx:8.3f}")

# Example 7: Signal modulation examples
print("\n=== Signal Modulation ===")

def amplitude_modulated_signal(t, carrier_freq=5, mod_freq=1, mod_depth=0.5):
    """AM signal: carrier modulated by triangle wave."""
    carrier = mp.sin(2 * mp.pi * carrier_freq * t)
    modulation = mp.trianglew(t, period=1/mod_freq)
    return (1 + mod_depth * modulation) * carrier

def pulse_width_modulated(t, freq=1, duty_cycle_func=None):
    """PWM signal with varying duty cycle."""
    if duty_cycle_func is None:
        duty_cycle_func = lambda x: 0.5 + 0.3 * mp.sin(2 * mp.pi * 0.2 * x)
    
    period = 1 / freq
    t_mod = t % period
    duty = duty_cycle_func(t)
    return 1 if t_mod < duty * period else -1

# Test modulated signals
t_mod = 0.5
am_signal = amplitude_modulated_signal(t_mod)
pwm_signal = pulse_width_modulated(t_mod)
print(f"AM signal at t={t_mod}: {am_signal:.4f}")
print(f"PWM signal at t={t_mod}: {pwm_signal}")

# Example 8: Window functions (using sigmoid)
print("\n=== Window Functions ===")

def sigmoid_window(t, width=2, steepness=5):
    """Smooth window using sigmoid functions."""
    left_edge = mp.sigmoid(steepness * (t + width/2))
    right_edge = mp.sigmoid(-steepness * (t - width/2))
    return left_edge * right_edge

window_t = [-2, -1, 0, 1, 2]
print("Time\tSigmoid Window")
for t in window_t:
    window_val = sigmoid_window(t, width=2, steepness=3)
    print(f"{t:4.1f}\t{window_val:.6f}")

# Example 9: Digital filter simulation
print("\n=== Digital Filter Simulation ===")

def moving_average_filter(signal_func, t, window_size=0.2, samples=21):
    """Simple moving average filter."""
    dt = window_size / (samples - 1)
    sum_val = 0
    for i in range(samples):
        sample_time = t - window_size/2 + i * dt
        sum_val += signal_func(sample_time)
    return sum_val / samples

# Apply filter to noisy square wave
def noisy_square(t):
    base = mp.squarew(t, period=1)
    noise = 0.1 * mp.sin(20 * mp.pi * t)  # High frequency noise
    return base + noise

test_times = [0, 0.1, 0.25, 0.4, 0.5]
print("Time\tNoisy Signal\tFiltered")
for t in test_times:
    noisy = noisy_square(t)
    filtered = moving_average_filter(noisy_square, t)
    print(f"{t:.2f}\t{noisy:8.4f}\t{filtered:8.4f}")

# Example 10: Harmonic analysis
print("\n=== Harmonic Analysis ===")

def analyze_harmonics(waveform_func, period=1, max_harmonic=5):
    """Simple harmonic analysis using integration."""
    fundamental_freq = 1 / period
    
    print(f"Harmonic analysis for period = {period}")
    print("n\tAmplitude\tPhase")
    
    for n in range(1, max_harmonic + 1):
        # Compute Fourier coefficients
        def integrand_cos(t):
            return waveform_func(t) * mp.cos(2 * mp.pi * n * fundamental_freq * t)
        
        def integrand_sin(t):
            return waveform_func(t) * mp.sin(2 * mp.pi * n * fundamental_freq * t)
        
        # Simple numerical integration (Riemann sum)
        a_n = 0
        b_n = 0
        samples = 100
        dt = period / samples
        
        for i in range(samples):
            t_sample = i * dt
            a_n += integrand_cos(t_sample) * dt
            b_n += integrand_sin(t_sample) * dt
        
        a_n *= 2 / period
        b_n *= 2 / period
        
        amplitude = mp.sqrt(a_n**2 + b_n**2)
        phase = mp.atan2(b_n, a_n) * 180 / mp.pi
        
        print(f"{n}\t{amplitude:.4f}\t\t{phase:.1f}°")

# Analyze square wave harmonics
analyze_harmonics(lambda t: mp.squarew(t, period=1), period=1, max_harmonic=3)

print("\n=== Signal Processing Examples Complete ===")