# -*- coding: utf-8 -*-
"""
Created on Sun Mar  9 19:39:13 2025

@author: AKourgli
"""

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import lfilter

def rcosdesign(beta, span, sps, norm='sqrt'):
    t = np.linspace(-span/2, span/2, span*sps)
    h = np.zeros_like(t)
    eps = 1e-9
    
    for i, ti in enumerate(t):
        if abs(ti) < 1e-8:
            h[i] = 1.0
        else:
            numerator = np.sin(np.pi*ti*(1 - beta)) + 4*beta*ti*np.cos(np.pi*ti*(1 + beta))
            denominator = np.pi*ti*(1 - (4*beta*ti)**2)
            
            if abs(denominator) < eps:
                h[i] = beta * np.sin(np.pi/(2*beta)) / 2
            else:
                h[i] = numerator / denominator
    
    h = np.maximum(h, 0)
    h /= np.linalg.norm(h)
    
    if norm == 'sqrt':
        h = np.sqrt(h)
    
    return h

# Paramètres de simulation
num_symbols = 200
sps = 8
beta = 0.25
additional_delay = 2
snr_db = 20  # Contrôle du niveau de bruit

# Génération des symboles avec préambule
symbols = np.random.randint(0, 2, num_symbols + 10)*2 - 1

# Création des signaux de référence
t = np.arange(len(symbols)*sps)/sps

# 1. Signal rectangulaire brut
rect_waveform = np.repeat(symbols, sps)

# 2. Filtrage RRC
up_samples = np.zeros_like(rect_waveform)
up_samples[::sps] = symbols
taps = rcosdesign(beta, 4, sps)
signal_shaped = lfilter(taps, 1.0, up_samples)

# 3. Application de l'ISI et calcul du retard
isi_taps = [0.5, 1.0, 0.5]
signal_isi = lfilter(isi_taps, 1.0, signal_shaped)

# Calcul des retards
group_delay_rrc = (len(taps)-1)//2
group_delay_isi = (len(isi_taps)-1)//2
total_delay = group_delay_rrc + group_delay_isi + additional_delay

# Ajout du bruit Gaussien
signal_power = np.mean(signal_isi**2)
noise_power = signal_power / (10 ** (snr_db / 10))
noise = np.random.normal(0, np.sqrt(noise_power), signal_isi.shape)
signal_noisy = signal_isi + noise


# Découpage du signal bruité
signal_noisy = signal_noisy[total_delay:]

# Découpage et alignement temporel
signal_shaped_cut = signal_shaped[total_delay:]  # Même découpe que le signal bruité
time_axis = np.arange(len(signal_shaped_cut))/sps  # Échelle temporelle correcte

# Préparation de l'affichage
plt.figure(figsize=(14, 10))

# 1. Symboles originaux
ax1 = plt.subplot(311)
ax1.stem(symbols, linefmt='C0-', markerfmt='C0o', basefmt=' ')
ax1.set_title('Symboles Transmis')
ax1.set_xlim(0, 50)
ax1.grid(True)

# 2. Comparaison des signaux
ax2 = plt.subplot(312)
ax2.plot(time_axis, signal_noisy, 'm', alpha=0.3, label=f'Signal bruité (SNR={snr_db}dB)')
ax2.plot(time_axis, signal_shaped_cut, 'g', alpha=0.8, label='Signal RRC original')
ax2.set_title('Comparaison des Signaux Avant/Après Bruit')
ax2.set_xlim(0, 200)
ax2.legend()
ax2.grid(True)

# 3. Diagramme de l'œil avec bruit
ax3 = plt.subplot(313)
eye_length = 2*sps
num_eyes = len(signal_noisy)//eye_length
t_eye = np.linspace(-0.5, 1.5, eye_length)

for i in range(num_eyes):
    segment = signal_noisy[i*eye_length:(i+1)*eye_length]
    ax3.plot(t_eye, segment, 'navy', alpha=0.5, linewidth=1)

ax3.set_title('Diagramme de l\'œil avec ISI et Bruit')
ax3.set_xlabel('Temps (Période Symbole)')
ax3.set_ylabel('Amplitude')
ax3.grid(True)
ax3.set_xlim(-0.5, 1.5)

plt.tight_layout()
plt.show()