Analiza wpływy usuniecia wyższych harmonicznych w zasilaniu silnika BLCD sygnałem PWM¶

przykłady powstały przy dużym wsparciu AI. Dokładnie Copilot-a GITHUB-a w środowsiku VSCode.

In [1]:
# Import bibliotek
import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots
import plotly.express as px

# musi być jeżeli chcemy do html-a
import plotly as plt
from IPython.display import display
plt.io.renderers.default = "notebook_connected" 

print("✅ Wykresy Plotly - w pełni interaktywne!")
print("  • Hover - najedź myszką aby zobaczyć wartości")
print("  • Zoom - zaznacz obszar lub użyj scroll")
print("  • Pan - przeciągnij wykres")
print("  • Kliknij legendę aby włączyć/wyłączyć serie")
print("  • Double-click aby zresetować widok")
print("  • Kliknij ikonę aparatu aby zapisać jako PNG")

✅ Wykresy Plotly - w pełni interaktywne!
  • Hover - najedź myszką aby zobaczyć wartości
  • Zoom - zaznacz obszar lub użyj scroll
  • Pan - przeciągnij wykres
  • Kliknij legendę aby włączyć/wyłączyć serie
  • Double-click aby zresetować widok
  • Kliknij ikonę aparatu aby zapisać jako PNG

Parametry podstawowe¶

In [2]:
# Parametry
V_max = 24         # Maksymalne napięcie fazowe [V]
f_electrical = 50  # Częstotliwość elektryczna [Hz] - symulowana prędkość
num_rotations = 2  # Liczba obrotów elektrycznych do pokazania

# OÅ› czasu
T_period = 1 / f_electrical  # Okres jednego obrotu [s]
t_end = num_rotations * T_period
dt = T_period / 1000  # 1000 próbek na obrót
t = np.arange(0, t_end, dt)

# KÄ…t elektryczny w czasie
theta_elec = 2 * np.pi * f_electrical * t  # [rad]

print(f"Parametry symulacji:")
print(f"  Napięcie maksymalne: {V_max} V")
print(f"  Częstotliwość: {f_electrical} Hz")
print(f"  Liczba obrotów: {num_rotations}")
print(f"  Czas symulacji: {t_end*1000:.1f} ms")
print(f"  Liczba próbek: {len(t)}")

Parametry symulacji:
  Napięcie maksymalne: 24 V
  Częstotliwość: 50 Hz
  Liczba obrotów: 2
  Czas symulacji: 40.0 ms
  Liczba próbek: 2000

Napięcia sinusoidalne na trzech fazach¶

Napięcia fazowe przesunięte o 120°:

  • Faza A: $V_A = V_{max} \cdot \sin(\theta)$
  • Faza B: $V_B = V_{max} \cdot \sin(\theta - \frac{2\pi}{3})$
  • Faza C: $V_C = V_{max} \cdot \sin(\theta - \frac{4\pi}{3})$
In [3]:
# Napięcia sinusoidalne na trzech fazach
V_A = V_max * np.sin(theta_elec)
V_B = V_max * np.sin(theta_elec - 2*np.pi/3)  # Przesunięcie o 120°
V_C = V_max * np.sin(theta_elec - 4*np.pi/3)  # Przesunięcie o 240°

# Weryfikacja - suma napięć powinna być zawsze 0
V_sum = V_A + V_B + V_C
print(f"Suma napięć (powinna być ~0): max = {np.max(np.abs(V_sum)):.2e} V")

Suma napięć (powinna być ~0): max = 2.49e-14 V

Wykres napięć w funkcji czasu¶

In [4]:
# Wykres napięć w czasie
fig = go.Figure()

fig.add_trace(go.Scatter(x=t * 1000, y=V_A, mode='lines', name='Faza A',
                         line=dict(color='red', width=2), opacity=0.8))
fig.add_trace(go.Scatter(x=t * 1000, y=V_B, mode='lines', name='Faza B',
                         line=dict(color='green', width=2), opacity=0.8))
fig.add_trace(go.Scatter(x=t * 1000, y=V_C, mode='lines', name='Faza C',
                         line=dict(color='blue', width=2), opacity=0.8))

# Linia zerowa
fig.add_hline(y=0, line_dash="solid", line_color="black", line_width=1, opacity=0.3)

# Zaznaczenie pełnych obrotów
for i in range(num_rotations + 1):
    t_mark = i * T_period * 1000
    fig.add_vline(x=t_mark, line_dash="dash", line_color="gray", line_width=1, opacity=0.4)
    if i > 0:
        fig.add_annotation(x=t_mark - T_period*500, y=V_max + 3, 
                          text=f'Obrót {i}', showarrow=False,
                          bgcolor="wheat", opacity=0.5)

fig.update_layout(
    title='Napięcia sinusoidalne na trzech fazach silnika BLDC',
    xaxis_title='Czas [ms]',
    yaxis_title='Napięcie [V]',
    hovermode='x unified',
    template='plotly_white',
    height=500,
    showlegend=True
)

fig.update_yaxes(range=[-V_max - 5, V_max + 5])

fig.show()

print(f"\nWłaściwości napięć:")
print(f"  Przesunięcie fazowe między fazami: 120° (2π/3 rad)")
print(f"  Amplituda: {V_max} V")
print(f"  Wartość RMS: {V_max/np.sqrt(2):.2f} V")

Właściwości napięć:
  Przesunięcie fazowe między fazami: 120° (2π/3 rad)
  Amplituda: 24 V
  Wartość RMS: 16.97 V

Wykres napięć w funkcji kąta obrotu¶

In [5]:
# Wykres napięć vs kąt
theta_deg = theta_elec * 180 / np.pi  # Konwersja na stopnie

fig = go.Figure()

fig.add_trace(go.Scatter(x=theta_deg, y=V_A, mode='lines', name='Faza A',
                         line=dict(color='red', width=2.5), opacity=0.8))
fig.add_trace(go.Scatter(x=theta_deg, y=V_B, mode='lines', name='Faza B',
                         line=dict(color='green', width=2.5), opacity=0.8))
fig.add_trace(go.Scatter(x=theta_deg, y=V_C, mode='lines', name='Faza C',
                         line=dict(color='blue', width=2.5), opacity=0.8))

# Linia zerowa
fig.add_hline(y=0, line_dash="solid", line_color="black", line_width=1, opacity=0.3)

# Zaznaczenie pełnych obrotów (360°)
for i in range(num_rotations + 1):
    angle_mark = i * 360
    fig.add_vline(x=angle_mark, line_dash="dash", line_color="gray", line_width=1.5, opacity=0.5)

# Zaznaczenie kluczowych punktów (co 60°)
for angle in range(0, num_rotations * 360 + 1, 60):
    fig.add_vline(x=angle, line_dash="dot", line_color="lightgray", line_width=0.8, opacity=0.4)

# Zaznaczenie punktów charakterystycznych
characteristic_angles = [0, 90, 180, 270, 360]
for angle in characteristic_angles:
    if angle <= 360:
        fig.add_annotation(x=angle, y=-V_max - 3, text=f'{angle}°', 
                          showarrow=False, font=dict(color='gray', size=9))

fig.update_layout(
    title='Napięcia fazowe vs. kąt obrotu rotora',
    xaxis_title='Kąt elektryczny [°]',
    yaxis_title='Napięcie [V]',
    hovermode='x unified',
    template='plotly_white',
    height=500,
    showlegend=True
)

fig.update_xaxes(range=[0, num_rotations * 360])
fig.update_yaxes(range=[-V_max - 5, V_max + 5])

fig.show()

print(f"\nObserwacje:")
print(f"  • Każda faza przechodzi przez pełny cykl sinusoidalny podczas 1 obrotu (360°)")
print(f"  • Fazy są przesunięte względem siebie o 120°")
print(f"  • W każdej chwili suma trzech napięć = 0 (układ symetryczny)")

Obserwacje:
  • Każda faza przechodzi przez pełny cykl sinusoidalny podczas 1 obrotu (360°)
  • Fazy są przesunięte względem siebie o 120°
  • W każdej chwili suma trzech napięć = 0 (układ symetryczny)

Wizualizacja 3D - wektory przestrzenne¶

In [6]:
# Reprezentacja fazorowa - pokaz kilka punktów z pierwszego obrotu
num_vectors = 12  # Pokaz 12 pozycji
angles_to_show = np.linspace(0, 2*np.pi, num_vectors, endpoint=False)

# Osie dla trzech faz (przesunięte o 120°)
phase_angles = [0, 2*np.pi/3, 4*np.pi/3]
phase_colors = ['red', 'green', 'blue']
phase_labels = ['Faza A (0°)', 'Faza B (120°)', 'Faza C (240°)']

fig = go.Figure()

# Narysuj osie faz
for angle, color, label in zip(phase_angles, phase_colors, phase_labels):
    fig.add_trace(go.Scatterpolar(
        r=[0, V_max], theta=[angle * 180/np.pi, angle * 180/np.pi],
        mode='lines', name=label, line=dict(color=color, width=3), opacity=0.3
    ))
    fig.add_trace(go.Scatterpolar(
        r=[0, V_max], theta=[(angle + np.pi) * 180/np.pi, (angle + np.pi) * 180/np.pi],
        mode='lines', name=label + ' (-)', line=dict(color=color, width=3, dash='dash'), 
        opacity=0.3, showlegend=False
    ))

# Pokaż wartości napięć w wybranych punktach
r_A, theta_A = [], []
r_B, theta_B = [], []
r_C, theta_C = [], []

for i, theta in enumerate(angles_to_show):
    V_A_val = V_max * np.sin(theta)
    V_B_val = V_max * np.sin(theta - 2*np.pi/3)
    V_C_val = V_max * np.sin(theta - 4*np.pi/3)
    
    # Faza A
    if V_A_val >= 0:
        r_A.append(abs(V_A_val))
        theta_A.append(phase_angles[0] * 180/np.pi)
    else:
        r_A.append(abs(V_A_val))
        theta_A.append((phase_angles[0] + np.pi) * 180/np.pi)
    
    # Faza B
    if V_B_val >= 0:
        r_B.append(abs(V_B_val))
        theta_B.append(phase_angles[1] * 180/np.pi)
    else:
        r_B.append(abs(V_B_val))
        theta_B.append((phase_angles[1] + np.pi) * 180/np.pi)
    
    # Faza C
    if V_C_val >= 0:
        r_C.append(abs(V_C_val))
        theta_C.append(phase_angles[2] * 180/np.pi)
    else:
        r_C.append(abs(V_C_val))
        theta_C.append((phase_angles[2] + np.pi) * 180/np.pi)

fig.add_trace(go.Scatterpolar(r=r_A, theta=theta_A, mode='markers',
                               marker=dict(color='red', size=10), opacity=0.6,
                               name='Próbki A', showlegend=False))
fig.add_trace(go.Scatterpolar(r=r_B, theta=theta_B, mode='markers',
                               marker=dict(color='green', size=10), opacity=0.6,
                               name='Próbki B', showlegend=False))
fig.add_trace(go.Scatterpolar(r=r_C, theta=theta_C, mode='markers',
                               marker=dict(color='blue', size=10), opacity=0.6,
                               name='Próbki C', showlegend=False))

fig.update_layout(
    title='Reprezentacja fazorowa napięć (widok biegunowy)<br>3 osie faz rozstawione co 120°',
    polar=dict(radialaxis=dict(range=[0, V_max + 2])),
    height=700,
    showlegend=True
)

fig.show()

print(f"\nDiagram biegunowy pokazuje:")
print(f"  • 3 osie faz rozstawione co 120°")
print(f"  • Napięcia dodatnie - na zewnątrz od środka")
print(f"  • Napięcia ujemne - w przeciwnym kierunku (poza 180°)")

Diagram biegunowy pokazuje:
  • 3 osie faz rozstawione co 120°
  • Napięcia dodatnie - na zewnątrz od środka
  • Napięcia ujemne - w przeciwnym kierunku (poza 180°)

Podsumowanie - Krok 1¶

✅ Co pokazaliśmy:

  • NapiÄ™cia na trzech fazach majÄ… ksztaÅ‚t sinusoidalny
  • SÄ… przesuniÄ™te wzglÄ™dem siebie o 120° (2Ï€/3 rad)
  • Suma trzech napięć w każdej chwili = 0 (ukÅ‚ad symetryczny)
  • Podczas 1 obrotu elektrycznego każda faza przechodzi przez peÅ‚ny cykl sin

Następny krok: Dodanie modulacji PWM do wygenerowania tych napięć sinusoidalnych.

Krok 2: Modulacja PWM¶

Sinusoidal PWM (SPWM)¶

Aby wygenerować napięcia sinusoidalne używając sterownika z tranzystorami (które mogą być tylko ON/OFF), stosujemy modulację PWM:

  • SygnaÅ‚ referencyjny: sinusoida (amplituda 0-1)
  • NoÅ›na PWM: sygnaÅ‚ trójkÄ…tny o wysokiej czÄ™stotliwoÅ›ci (20 kHz)
  • WyjÅ›cie: PWM = 1 gdy sinus > noÅ›na, PWM = 0 w przeciwnym przypadku

Parametry PWM¶

In [7]:
# Parametry PWM
f_pwm = 20000      # Częstotliwość PWM [Hz]
modulation_index = 0.9  # Indeks modulacji (0-1) - kontroluje amplitudÄ™

# Nowy wektor czasu z dużo wyższą rozdzielczością dla PWM
dt_pwm = 1 / (f_pwm * 100)  # 100 próbek na okres PWM
t_pwm = np.arange(0, t_end, dt_pwm)

# KÄ…t elektryczny dla nowego wektora czasu
theta_elec_pwm = 2 * np.pi * f_electrical * t_pwm

# Napięcia sinusoidalne znormalizowane (0 do 1)
V_A_norm = modulation_index * np.sin(theta_elec_pwm)
V_B_norm = modulation_index * np.sin(theta_elec_pwm - 2*np.pi/3)
V_C_norm = modulation_index * np.sin(theta_elec_pwm - 4*np.pi/3)

# Sygnał nośnej PWM (trójkątny od -1 do +1)
T_pwm = 1 / f_pwm
carrier = 2 * np.abs(2 * ((t_pwm % T_pwm) / T_pwm) - 1) - 1  # Trójkątny -1 do +1

print(f"Parametry modulacji PWM:")
print(f"  Częstotliwość PWM: {f_pwm/1000} kHz")
print(f"  Indeks modulacji: {modulation_index}")
print(f"  Liczba próbek: {len(t_pwm)}")
print(f"  Rozdzielczość czasowa: {dt_pwm*1e6:.2f} μs")
print(f"  Stosunek częstotliwości PWM/elektrycznej: {f_pwm/f_electrical:.0f}")

Parametry modulacji PWM:
  Częstotliwość PWM: 20.0 kHz
  Indeks modulacji: 0.9
  Liczba próbek: 80000
  Rozdzielczość czasowa: 0.50 μs
  Stosunek częstotliwości PWM/elektrycznej: 400

Generowanie sygnałów PWM¶

In [8]:
# Napięcie zasilania DC
V_dc = V_max  # Napięcie DC [V]

# Komparator: PWM = 1 gdy sinus > carrier, PWM = 0 w przeciwnym razie
PWM_A = (V_A_norm > carrier).astype(float)
PWM_B = (V_B_norm > carrier).astype(float)
PWM_C = (V_C_norm > carrier).astype(float)

# Napięcia wyjściowe (z pełnym zakresem -V_dc do +V_dc)
# PWM = 0 → -V_dc, PWM = 1 → +V_dc
V_A_pwm = V_dc * (2 * PWM_A - 1)
V_B_pwm = V_dc * (2 * PWM_B - 1)
V_C_pwm = V_dc * (2 * PWM_C - 1)

# Średnie wartości (po filtracji)
V_A_target = V_dc * modulation_index * np.sin(theta_elec_pwm)
V_B_target = V_dc * modulation_index * np.sin(theta_elec_pwm - 2*np.pi/3)
V_C_target = V_dc * modulation_index * np.sin(theta_elec_pwm - 4*np.pi/3)

print(f"Sygnały PWM wygenerowane dla trzech faz")
print(f"  Średnie wypełnienie PWM dla fazy A: {np.mean(PWM_A)*100:.1f}%")
print(f"  Średnie wypełnienie PWM dla fazy B: {np.mean(PWM_B)*100:.1f}%")
print(f"  Średnie wypełnienie PWM dla fazy C: {np.mean(PWM_C)*100:.1f}%")

Sygnały PWM wygenerowane dla trzech faz
  Średnie wypełnienie PWM dla fazy A: 50.0%
  Średnie wypełnienie PWM dla fazy B: 50.0%
  Średnie wypełnienie PWM dla fazy C: 50.0%

Wizualizacja zasady działania SPWM (zoom na fragment)¶

In [9]:
# Zoom na krótki fragment czasu aby pokazać zasadę PWM
t_zoom_start = 0.005  # 5 ms
t_zoom_end = 0.008    # 8 ms (3 ms okno)
zoom_idx_start = np.searchsorted(t_pwm, t_zoom_start)
zoom_idx_end = np.searchsorted(t_pwm, t_zoom_end)

t_zoom = t_pwm[zoom_idx_start:zoom_idx_end]
V_A_norm_zoom = V_A_norm[zoom_idx_start:zoom_idx_end]
carrier_zoom = carrier[zoom_idx_start:zoom_idx_end]
PWM_A_zoom = PWM_A[zoom_idx_start:zoom_idx_end]
V_A_pwm_zoom = V_A_pwm[zoom_idx_start:zoom_idx_end]

# Wykres zasady działania PWM
fig = make_subplots(rows=2, cols=1, 
                    subplot_titles=('Zasada modulacji SPWM - Faza A (sygnał referencyjny vs nośna)',
                                   'Wyjściowe napięcie PWM na fazie A'),
                    vertical_spacing=0.12)

# Porównanie sinusa z nośną
fig.add_trace(go.Scatter(x=t_zoom * 1000, y=V_A_norm_zoom, mode='lines',
                         name='Sygnał referencyjny (sinus)',
                         line=dict(color='red', width=2.5), opacity=0.9),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_zoom * 1000, y=carrier_zoom, mode='lines',
                         name='Nośna PWM (trójkątna, 20 kHz)',
                         line=dict(color='blue', width=1.5), opacity=0.7),
              row=1, col=1)

# Wypełnienie obszarów
# PWM ON
mask_on = V_A_norm_zoom >= carrier_zoom
if np.any(mask_on):
    fig.add_trace(go.Scatter(x=t_zoom[mask_on] * 1000, y=V_A_norm_zoom[mask_on],
                            fill='tonexty', mode='none', name='PWM = 1 (ON)',
                            fillcolor='rgba(0,255,0,0.3)', showlegend=True),
                 row=1, col=1)

# Wyjściowy sygnał PWM
fig.add_trace(go.Scatter(x=t_zoom * 1000, y=V_A_pwm_zoom, mode='lines',
                         name='Napięcie PWM na fazie A',
                         line=dict(color='purple', width=1), opacity=0.8),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_zoom * 1000, y=V_A_target[zoom_idx_start:zoom_idx_end],
                         mode='lines', name='Wartość średnia (cel)',
                         line=dict(color='red', width=2.5, dash='dash'), opacity=0.9),
              row=2, col=1)

fig.update_xaxes(title_text='Czas [ms]', row=1, col=1)
fig.update_xaxes(title_text='Czas [ms]', row=2, col=1)
fig.update_yaxes(title_text='Amplituda znormalizowana', range=[-1.2, 1.2], row=1, col=1)
fig.update_yaxes(title_text='Napięcie [V]', range=[-V_dc - 5, V_dc + 5], row=2, col=1)

fig.update_layout(height=700, hovermode='x unified', template='plotly_white')

fig.show()

print(f"\nZasada działania:")
print(f"  • Gdy sinus > nośna → PWM = ON → napięcie = +{V_dc}V")
print(f"  • Gdy sinus < nośna → PWM = OFF → napięcie = -{V_dc}V")
print(f"  • Średnia wartość PWM = wartość sinusa × {V_dc}V")

Zasada działania:
  • Gdy sinus > nośna → PWM = ON → napięcie = +24V
  • Gdy sinus < nośna → PWM = OFF → napięcie = -24V
  • Średnia wartość PWM = wartość sinusa × 24V

Napięcia PWM na trzech fazach - cały przebieg¶

In [10]:
# Wykres napięć PWM na wszystkich fazach (pełny przebieg)
fig = make_subplots(rows=2, cols=1,
                    subplot_titles=('Napięcia PWM na trzech fazach - pełny przebieg',
                                   'Wartości średnie napięć PWM (po filtracji przez indukcyjność silnika)'),
                    vertical_spacing=0.12)

# Napięcia PWM (zbyt gęste dla pełnego wykresu - pokaz co 10-tą próbkę)
decimation = 10
t_plot = t_pwm[::decimation]
V_A_pwm_plot = V_A_pwm[::decimation]
V_B_pwm_plot = V_B_pwm[::decimation]
V_C_pwm_plot = V_C_pwm[::decimation]

fig.add_trace(go.Scatter(x=t_plot * 1000, y=V_A_pwm_plot, mode='lines',
                         name='Faza A', line=dict(color='red', width=0.3), opacity=0.7),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_plot * 1000, y=V_B_pwm_plot, mode='lines',
                         name='Faza B', line=dict(color='green', width=0.3), opacity=0.7),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_plot * 1000, y=V_C_pwm_plot, mode='lines',
                         name='Faza C', line=dict(color='blue', width=0.3), opacity=0.7),
              row=1, col=1)

# Wartości średnie (docelowe sinusoidy)
fig.add_trace(go.Scatter(x=t_pwm * 1000, y=V_A_target, mode='lines',
                         name='Faza A (średnia)', line=dict(color='red', width=2), opacity=0.8,
                         showlegend=False),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_pwm * 1000, y=V_B_target, mode='lines',
                         name='Faza B (średnia)', line=dict(color='green', width=2), opacity=0.8,
                         showlegend=False),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_pwm * 1000, y=V_C_target, mode='lines',
                         name='Faza C (średnia)', line=dict(color='blue', width=2), opacity=0.8,
                         showlegend=False),
              row=2, col=1)

# Zaznaczenie pełnych obrotów
for i in range(num_rotations + 1):
    t_mark = i * T_period * 1000
    fig.add_vline(x=t_mark, line_dash="dash", line_color="gray", line_width=1, opacity=0.4, row=2, col=1)

# Linie zerowe
fig.add_hline(y=0, line_dash="solid", line_color="black", line_width=1, opacity=0.3, row=1, col=1)
fig.add_hline(y=0, line_dash="solid", line_color="black", line_width=1, opacity=0.3, row=2, col=1)

fig.update_xaxes(title_text='Czas [ms]', row=1, col=1)
fig.update_xaxes(title_text='Czas [ms]', row=2, col=1)
fig.update_yaxes(title_text='Napięcie [V]', range=[-V_dc - 5, V_dc + 5], row=1, col=1)
fig.update_yaxes(title_text='Napięcie [V]', range=[-V_dc - 5, V_dc + 5], row=2, col=1)

fig.update_layout(height=700, hovermode='x unified', template='plotly_white')

fig.show()

print(f"Napięcia PWM:")
print(f"  Maksymalne napięcie: ±{V_dc} V")
print(f"  Napięcie średnie maksymalne: ±{V_dc * modulation_index:.1f} V")
print(f"  Częstotliwość przełączeń: {f_pwm/1000} kHz")

Napięcia PWM:
  Maksymalne napięcie: ±24 V
  Napięcie średnie maksymalne: ±21.6 V
  Częstotliwość przełączeń: 20.0 kHz

Zoom na jeden obrót - widok PWM¶

In [11]:
# Zoom na fragment jednego obrotu
t_zoom2_start = T_period * 0.5  # Åšrodek pierwszego obrotu
t_zoom2_end = t_zoom2_start + 0.003  # 3 ms
zoom2_idx_start = np.searchsorted(t_pwm, t_zoom2_start)
zoom2_idx_end = np.searchsorted(t_pwm, t_zoom2_end)

fig = go.Figure()

# Wszystkie trzy fazy - PWM
fig.add_trace(go.Scatter(x=t_pwm[zoom2_idx_start:zoom2_idx_end] * 1000,
                         y=V_A_pwm[zoom2_idx_start:zoom2_idx_end],
                         mode='lines', name='Faza A (PWM)',
                         line=dict(color='red', width=0.8), opacity=0.7))
fig.add_trace(go.Scatter(x=t_pwm[zoom2_idx_start:zoom2_idx_end] * 1000,
                         y=V_B_pwm[zoom2_idx_start:zoom2_idx_end],
                         mode='lines', name='Faza B (PWM)',
                         line=dict(color='green', width=0.8), opacity=0.7))
fig.add_trace(go.Scatter(x=t_pwm[zoom2_idx_start:zoom2_idx_end] * 1000,
                         y=V_C_pwm[zoom2_idx_start:zoom2_idx_end],
                         mode='lines', name='Faza C (PWM)',
                         line=dict(color='blue', width=0.8), opacity=0.7))

# Wartości średnie
fig.add_trace(go.Scatter(x=t_pwm[zoom2_idx_start:zoom2_idx_end] * 1000,
                         y=V_A_target[zoom2_idx_start:zoom2_idx_end],
                         mode='lines', name='Faza A (średnia)',
                         line=dict(color='red', width=2.5, dash='dash'), opacity=0.9))
fig.add_trace(go.Scatter(x=t_pwm[zoom2_idx_start:zoom2_idx_end] * 1000,
                         y=V_B_target[zoom2_idx_start:zoom2_idx_end],
                         mode='lines', name='Faza B (średnia)',
                         line=dict(color='green', width=2.5, dash='dash'), opacity=0.9))
fig.add_trace(go.Scatter(x=t_pwm[zoom2_idx_start:zoom2_idx_end] * 1000,
                         y=V_C_target[zoom2_idx_start:zoom2_idx_end],
                         mode='lines', name='Faza C (średnia)',
                         line=dict(color='blue', width=2.5, dash='dash'), opacity=0.9))

fig.add_hline(y=0, line_dash="solid", line_color="black", line_width=1, opacity=0.3)

fig.update_layout(
    title='Napięcia PWM trzech faz - zoom (widoczna modulacja 20 kHz)',
    xaxis_title='Czas [ms]',
    yaxis_title='Napięcie [V]',
    hovermode='x unified',
    template='plotly_white',
    height=500,
    showlegend=True
)

fig.update_yaxes(range=[-V_dc - 5, V_dc + 5])

fig.show()

print(f"Na wykresie widać:")
print(f"  • Szybkie przełączenia PWM (20 kHz)")
print(f"  • Zmiana wypełnienia PWM zgodnie z sinusoidą")
print(f"  • Fazy przesunięte o 120°")

Na wykresie widać:
  • Szybkie przełączenia PWM (20 kHz)
  • Zmiana wypełnienia PWM zgodnie z sinusoidą
  • Fazy przesunięte o 120°

Podsumowanie - Krok 2¶

✅ Co dodaliśmy:

  • Modulacja PWM o czÄ™stotliwoÅ›ci 20 kHz
  • Porównanie sinusoidy z noÅ›nÄ… trójkÄ…tnÄ… (komparator)
  • WyjÅ›ciowe sygnaÅ‚y PWM na trzech fazach
  • WartoÅ›ci Å›rednie odpowiadajÄ… sinusoidom

Kluczowe:

  • PWM przeÅ‚Ä…cza siÄ™ miÄ™dzy +24V i -24V
  • WypeÅ‚nienie PWM zmienia siÄ™ sinusoidalnie
  • Åšrednia wartość (po filtracji przez indukcyjność) = sinus × 24V × 0.9
  • Indukcyjność uzwojeÅ„ silnika dziaÅ‚a jak filtr dolnoprzepustowy

Następny krok: Symulacja prądów w uzwojeniach silnika z indukcyjnością (filtracja PWM)

Krok 3: Analiza częstotliwościowa (FFT)¶

Widmo sygnałów PWM¶

Analiza FFT pokazuje składowe częstotliwościowe w sygnale PWM:

  • CzÄ™stotliwość podstawowa (50 Hz) - sygnaÅ‚ użyteczny
  • CzÄ™stotliwość PWM (20 kHz) i jej harmoniczne - szum przeÅ‚Ä…czeÅ„
  • Produkty intermodulacji
In [12]:
# FFT sygnału PWM - Faza A
# Używamy tylko części w stanie ustalonym (pomijamy początek)
start_idx = len(V_A_pwm) // 4  # Rozpocznij od 25%
V_A_pwm_fft = V_A_pwm[start_idx:]
t_pwm_fft = t_pwm[start_idx:] - t_pwm[start_idx]

# FFT
N = len(V_A_pwm_fft)
fft_values = np.fft.fft(V_A_pwm_fft)
fft_freq = np.fft.fftfreq(N, dt_pwm)

# Tylko dodatnie częstotliwości
positive_freq_idx = fft_freq > 0
fft_freq_positive = fft_freq[positive_freq_idx]
fft_magnitude = np.abs(fft_values[positive_freq_idx]) * 2 / N  # Amplituda

print(f"Analiza FFT:")
print(f"  Liczba próbek: {N}")
print(f"  Rozdzielczość częstotliwościowa: {1/(t_pwm_fft[-1]):.2f} Hz")
print(f"  Maksymalna częstotliwość (Nyquist): {1/(2*dt_pwm)/1000:.1f} kHz")

Analiza FFT:
  Liczba próbek: 60000
  Rozdzielczość częstotliwościowa: 33.33 Hz
  Maksymalna częstotliwość (Nyquist): 1000.0 kHz

Widmo częstotliwościowe - pełny zakres¶

In [13]:
# Wykres widma - skala logarytmiczna (pełny zakres)
fig = go.Figure()

fig.add_trace(go.Scatter(x=fft_freq_positive / 1000, y=fft_magnitude,
                         mode='lines', name='Widmo FFT',
                         line=dict(color='blue', width=1)))

# Zaznaczenie kluczowych częstotliwości
fig.add_vline(x=f_electrical / 1000, line_dash="dash", line_color="green", 
              line_width=2, annotation_text=f"f₀ = {f_electrical} Hz (sygnał użyteczny)")
fig.add_vline(x=f_pwm / 1000, line_dash="dash", line_color="red", 
              line_width=2, annotation_text=f"f_PWM = {f_pwm/1000} kHz")

# Harmoniczne PWM (2f, 3f, 4f)
for n in [2, 3, 4]:
    fig.add_vline(x=n * f_pwm / 1000, line_dash="dot", line_color="red", 
                  line_width=1, opacity=0.5, annotation_text=f"{n}×f_PWM")

fig.update_layout(
    title='Widmo częstotliwościowe sygnału PWM - Faza A (skala logarytmiczna)',
    xaxis_title='Częstotliwość [kHz]',
    yaxis_title='Amplituda [V]',
    yaxis_type='log',
    template='plotly_white',
    height=600,
    showlegend=True,
    hovermode='x unified'
)

fig.update_xaxes(range=[0, 100])  # Do 100 kHz

fig.show()

print(f"\nObserwacje:")
print(f"  • Szczyt przy {f_electrical} Hz - częstotliwość podstawowa (sygnał użyteczny)")
print(f"  • Szczyt przy {f_pwm/1000} kHz - częstotliwość PWM")
print(f"  • Widoczne harmoniczne PWM przy wielokrotnościach 20 kHz")

Obserwacje:
  • Szczyt przy 50 Hz - częstotliwość podstawowa (sygnał użyteczny)
  • Szczyt przy 20.0 kHz - częstotliwość PWM
  • Widoczne harmoniczne PWM przy wielokrotnościach 20 kHz

Zoom na niskie częstotliwości (sygnał użyteczny)¶

In [14]:
# Zoom na niskie częstotliwości (0-500 Hz)
freq_limit_low = 500  # Hz
low_freq_idx = fft_freq_positive < freq_limit_low
fft_freq_low = fft_freq_positive[low_freq_idx]
fft_magnitude_low = fft_magnitude[low_freq_idx]

fig = go.Figure()

fig.add_trace(go.Scatter(x=fft_freq_low, y=fft_magnitude_low,
                         mode='lines+markers', name='Widmo FFT (niskie częstotliwości)',
                         line=dict(color='blue', width=2),
                         marker=dict(size=4)))

# Zaznaczenie częstotliwości podstawowej
fig.add_vline(x=f_electrical, line_dash="dash", line_color="green", 
              line_width=2, annotation_text=f"fâ‚€ = {f_electrical} Hz")

# Harmoniczne częstotliwości podstawowej (3f, 5f, 7f - nieparzyste)
for n in [3, 5, 7]:
    if n * f_electrical < freq_limit_low:
        fig.add_vline(x=n * f_electrical, line_dash="dot", line_color="orange", 
                      line_width=1, opacity=0.5, annotation_text=f"{n}×f₀")

fig.update_layout(
    title='Widmo częstotliwościowe - Zoom na sygnał użyteczny (0-500 Hz)',
    xaxis_title='Częstotliwość [Hz]',
    yaxis_title='Amplituda [V]',
    template='plotly_white',
    height=600,
    showlegend=True,
    hovermode='x unified'
)

fig.show()

# Znajdź dominujące składowe
top_n = 10
top_indices = np.argsort(fft_magnitude_low)[-top_n:][::-1]
print(f"\nTop {top_n} składowych częstotliwościowych (0-{freq_limit_low} Hz):")
for idx in top_indices:
    freq = fft_freq_low[idx]
    amp = fft_magnitude_low[idx]
    if amp > 0.01:  # Próg
        print(f"  f = {freq:.1f} Hz, Amplituda = {amp:.2f} V")

Top 10 składowych częstotliwościowych (0-500 Hz):
  f = 33.3 Hz, Amplituda = 16.52 V
  f = 66.7 Hz, Amplituda = 11.80 V
  f = 100.0 Hz, Amplituda = 3.04 V
  f = 133.3 Hz, Amplituda = 1.50 V
  f = 166.7 Hz, Amplituda = 0.90 V
  f = 200.0 Hz, Amplituda = 0.61 V
  f = 233.3 Hz, Amplituda = 0.42 V
  f = 266.7 Hz, Amplituda = 0.35 V
  f = 300.0 Hz, Amplituda = 0.28 V
  f = 333.3 Hz, Amplituda = 0.22 V

Zoom na częstotliwość PWM i harmoniczne¶

In [15]:
# Zoom na pasmo PWM z harmonicznymi (0-100 kHz)
freq_low_pwm = 0  # 0 kHz
freq_high_pwm = 100000  # 100 kHz - obejmuje pierwsze 5 harmonicznych PWM
pwm_band_idx = (fft_freq_positive >= freq_low_pwm) & (fft_freq_positive <= freq_high_pwm)
fft_freq_pwm = fft_freq_positive[pwm_band_idx]
fft_magnitude_pwm = fft_magnitude[pwm_band_idx]

fig = go.Figure()

fig.add_trace(go.Scatter(x=fft_freq_pwm / 1000, y=fft_magnitude_pwm,
                         mode='lines', name='Widmo FFT (pasmo PWM z harmonicznymi)',
                         line=dict(color='red', width=1.5)))

# Zaznaczenie częstotliwości podstawowej
fig.add_vline(x=f_electrical / 1000, line_dash="dash", line_color="green", 
              line_width=2, annotation_text=f"fâ‚€ = {f_electrical} Hz", 
              annotation_position="top left")

# Harmoniczne PWM (1f, 2f, 3f, 4f, 5f)
colors_harmonics = ['red', 'orange', 'purple', 'brown', 'pink']
for n in range(1, 6):
    color = colors_harmonics[n-1]
    fig.add_vline(x=n * f_pwm / 1000, line_dash="dash", line_color=color, 
                  line_width=2, opacity=0.7, 
                  annotation_text=f"{n}×f_PWM = {n*f_pwm/1000} kHz",
                  annotation_position="top")

# Pasma boczne wokół pierwszej harmonicznej PWM (f_PWM ± f₀)
fig.add_vline(x=(f_pwm - f_electrical) / 1000, line_dash="dot", line_color="orange", 
              line_width=1, opacity=0.5)
fig.add_vline(x=(f_pwm + f_electrical) / 1000, line_dash="dot", line_color="orange", 
              line_width=1, opacity=0.5)

# Pasma boczne wokół drugiej harmonicznej PWM (2×f_PWM ± f₀)
fig.add_vline(x=(2*f_pwm - f_electrical) / 1000, line_dash="dot", line_color="orange", 
              line_width=1, opacity=0.4)
fig.add_vline(x=(2*f_pwm + f_electrical) / 1000, line_dash="dot", line_color="orange", 
              line_width=1, opacity=0.4)

# Pasma boczne wokół trzeciej harmonicznej PWM (3×f_PWM ± f₀)
fig.add_vline(x=(3*f_pwm - f_electrical) / 1000, line_dash="dot", line_color="orange", 
              line_width=1, opacity=0.3)
fig.add_vline(x=(3*f_pwm + f_electrical) / 1000, line_dash="dot", line_color="orange", 
              line_width=1, opacity=0.3)

fig.update_layout(
    title='Widmo częstotliwościowe - PWM z harmonicznymi (0-100 kHz)',
    xaxis_title='Częstotliwość [kHz]',
    yaxis_title='Amplituda [V]',
    yaxis_type='log',  # Skala logarytmiczna dla lepszej widoczności
    template='plotly_white',
    height=700,
    showlegend=True,
    hovermode='x unified'
)

fig.show()

print(f"\nW szerokim paśmie PWM widać:")
print(f"  • Częstotliwość podstawowa: {f_electrical} Hz")
print(f"  • 1. harmoniczna PWM: {f_pwm/1000} kHz (najsilniejsza)")
print(f"  • 2. harmoniczna PWM: {2*f_pwm/1000} kHz")
print(f"  • 3. harmoniczna PWM: {3*f_pwm/1000} kHz")
print(f"  • 4. harmoniczna PWM: {4*f_pwm/1000} kHz")
print(f"  • 5. harmoniczna PWM: {5*f_pwm/1000} kHz")
print(f"  • Pasma boczne (sidebands) wokół każdej harmonicznej przy n×f_PWM ± f₀")
print(f"  • Amplitudy harmonicznych maleją z rzędem (filtracja przez indukcyjność)")

W szerokim paśmie PWM widać:
  • Częstotliwość podstawowa: 50 Hz
  • 1. harmoniczna PWM: 20.0 kHz (najsilniejsza)
  • 2. harmoniczna PWM: 40.0 kHz
  • 3. harmoniczna PWM: 60.0 kHz
  • 4. harmoniczna PWM: 80.0 kHz
  • 5. harmoniczna PWM: 100.0 kHz
  • Pasma boczne (sidebands) wokół każdej harmonicznej przy n×f_PWM ± f₀
  • Amplitudy harmonicznych maleją z rzędem (filtracja przez indukcyjność)

Porównanie widm trzech faz¶

In [16]:
# FFT dla wszystkich trzech faz
V_B_pwm_fft = V_B_pwm[start_idx:]
V_C_pwm_fft = V_C_pwm[start_idx:]

fft_B = np.fft.fft(V_B_pwm_fft)
fft_C = np.fft.fft(V_C_pwm_fft)

fft_magnitude_B = np.abs(fft_B[positive_freq_idx]) * 2 / N
fft_magnitude_C = np.abs(fft_C[positive_freq_idx]) * 2 / N

# Porównanie (tylko niskie częstotliwości)
fig = go.Figure()

fig.add_trace(go.Scatter(x=fft_freq_low, y=fft_magnitude_low,
                         mode='lines', name='Faza A',
                         line=dict(color='red', width=2), opacity=0.8))
fig.add_trace(go.Scatter(x=fft_freq_low, y=fft_magnitude_B[low_freq_idx],
                         mode='lines', name='Faza B',
                         line=dict(color='green', width=2), opacity=0.8))
fig.add_trace(go.Scatter(x=fft_freq_low, y=fft_magnitude_C[low_freq_idx],
                         mode='lines', name='Faza C',
                         line=dict(color='blue', width=2), opacity=0.8))

fig.update_layout(
    title='Porównanie widm trzech faz (0-500 Hz)',
    xaxis_title='Częstotliwość [Hz]',
    yaxis_title='Amplituda [V]',
    template='plotly_white',
    height=600,
    showlegend=True,
    hovermode='x unified'
)

fig.show()

print(f"\nPorównanie faz:")
print(f"  • Wszystkie fazy mają identyczne widmo")
print(f"  • Różnią się tylko przesunięciem fazowym (120°)")
print(f"  • Amplitudy składowych są takie same")

Porównanie faz:
  • Wszystkie fazy mają identyczne widmo
  • Różnią się tylko przesunięciem fazowym (120°)
  • Amplitudy składowych są takie same

Podsumowanie - Analiza FFT¶

✅ Co pokazała analiza widmowa:

  1. Częstotliwość podstawowa (50 Hz):

    • NajwiÄ™ksza skÅ‚adowa - to nasz sygnaÅ‚ użyteczny
    • Tworzy pole magnetyczne obrotowe w silniku

  2. Częstotliwość PWM (20 kHz):

    • Widoczna jako dominujÄ…ca skÅ‚adowa wysokoczÄ™stotliwoÅ›ciowa
    • Wraz z harmonicznymi (40 kHz, 60 kHz...)
    • Te skÅ‚adowe sÄ… filtowane przez indukcyjność silnika

  3. Pasma boczne:

    • Wokół czÄ™stotliwoÅ›ci PWM (f_PWM ± fâ‚€)
    • Wynik intermodulacji miÄ™dzy sygnaÅ‚em a noÅ›nÄ…

  4. Znaczenie dla silnika:

    • Indukcyjność uzwojeÅ„ (L) dziaÅ‚a jak filtr dolnoprzepustowy
    • TÅ‚umi skÅ‚adowe o czÄ™stotliwoÅ›ci PWM
    • PrÄ…d ma głównie skÅ‚adowÄ… 50 Hz (gÅ‚adki, sinusoidalny)
    • MaÅ‚e tÄ™tnienia prÄ…du przy 20 kHz

Następny krok: Symulacja prądów w uzwojeniach z uwzględnieniem filtracji przez indukcyjność

In [17]:
# Skuteczny sinus z PWM: prosty filtr dolnoprzepustowy (1. rzÄ…d)
fc_lp = 1000  # [Hz] odcięcie: dużo wyżej niż 50 Hz, dużo niżej niż 20 kHz
RC = 1 / (2 * np.pi * fc_lp)
alpha = dt_pwm / (RC + dt_pwm)

V_A_filtered = np.zeros_like(V_A_pwm)
for i in range(1, len(V_A_pwm)):
    V_A_filtered[i] = V_A_filtered[i-1] + alpha * (V_A_pwm[i] - V_A_filtered[i-1])

# RMS i błąd względem sinusoidy referencyjnej
V_A_rms = np.sqrt(np.mean(V_A_filtered**2))
V_A_target_rms = np.sqrt(np.mean(V_A_target**2))
max_err = np.max(np.abs(V_A_filtered - V_A_target))

# Wizualizacja: pełny przebieg + zoom
decimation_lp = 20
t_plot_lp = t_pwm[::decimation_lp] * 1000

t_zoom_start_lp = 0.005
t_zoom_end_lp = 0.008
zoom_mask = (t_pwm >= t_zoom_start_lp) & (t_pwm <= t_zoom_end_lp)

fig = make_subplots(rows=2, cols=1, vertical_spacing=0.12,
                    subplot_titles=('PWM vs filtr LP (pełny przebieg)',
                                    'Zoom: PWM → sinus po filtracji'))

fig.add_trace(go.Scatter(x=t_plot_lp, y=V_A_pwm[::decimation_lp], mode='lines',
                        name='PWM (A)', line=dict(color='gray', width=0.6), opacity=0.6),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_plot_lp, y=V_A_filtered[::decimation_lp], mode='lines',
                        name='Po filtrze (A)', line=dict(color='red', width=2)),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_plot_lp, y=V_A_target[::decimation_lp], mode='lines',
                        name='Docelowy sinus (A)', line=dict(color='black', width=2, dash='dash')),
              row=1, col=1)

fig.add_trace(go.Scatter(x=t_pwm[zoom_mask] * 1000, y=V_A_pwm[zoom_mask], mode='lines',
                        name='PWM (zoom)', line=dict(color='gray', width=0.6), opacity=0.6),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_pwm[zoom_mask] * 1000, y=V_A_filtered[zoom_mask], mode='lines',
                        name='Po filtrze (zoom)', line=dict(color='red', width=2)),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_pwm[zoom_mask] * 1000, y=V_A_target[zoom_mask], mode='lines',
                        name='Docelowy sinus (zoom)', line=dict(color='black', width=2, dash='dash')),
              row=2, col=1)

fig.update_xaxes(title_text='Czas [ms]', row=1, col=1)
fig.update_xaxes(title_text='Czas [ms]', row=2, col=1)
fig.update_yaxes(title_text='Napięcie [V]', row=1, col=1)
fig.update_yaxes(title_text='Napięcie [V]', row=2, col=1)

fig.update_layout(template='plotly_white', height=700, hovermode='x unified')

fig.show()

print(f"V_A skuteczne po filtracji: {V_A_rms:.2f} V RMS")
print(f"Docelowy sinus RMS:         {V_A_target_rms:.2f} V RMS")
print(f"Maks. odchyłka chwilowa:    {max_err:.2f} V")

V_A skuteczne po filtracji: 15.29 V RMS
Docelowy sinus RMS:         15.27 V RMS
Maks. odchyłka chwilowa:    3.01 V

Powyżej widać, że AI podsunęła pomysł zastosowania filtra LP, aby zasymulować wygładzenie PWM na uzwojeniu. Same uzwojenie silnika działa jak filtr LR (cewka).¶

Sugestia ta potwierdza, że zastosowanie tego filtra jest najkorzystniejsze¶

Jeszcze inny typy filtrów¶

próba usuniecia harmonicznych filtrem notch¶

In [18]:
# Filtracja harmonicznych PWM metodą FFT (notch na wielokrotnościach f_pwm)
import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots

# Parametry filtra
stop_bw = 2000       # szerokość wycięcia wokół każdej harmonicznej [Hz]
max_harm = 6         # liczba wycinanych harmonicznych f_pwm (1..max_harm)

# FFT sygnału PWM fazy A
N = len(V_A_pwm)
freqs = np.fft.fftfreq(N, dt_pwm)
fft_vals = np.fft.fft(V_A_pwm)

# Maska: przepuszczamy wszystko poza pasmami n*f_pwm ± stop_bw
mask = np.ones_like(fft_vals, dtype=bool)
for n in range(1, max_harm + 1):
    band = np.abs(np.abs(freqs) - n * f_pwm) < stop_bw
    mask &= ~band

fft_filtered = fft_vals * mask
V_A_notch = np.fft.ifft(fft_filtered).real

# Metryki
rms_notch = np.sqrt(np.mean(V_A_notch**2))
rms_target = np.sqrt(np.mean(V_A_target**2))
max_err = np.max(np.abs(V_A_notch - V_A_target))

# Wykresy: pełny przebieg + zoom
zoom_start = 0.005
zoom_end = 0.008
zoom_mask = (t_pwm >= zoom_start) & (t_pwm <= zoom_end)

fig = make_subplots(rows=2, cols=1, vertical_spacing=0.12,
                    subplot_titles=(
                        "PWM vs filtr notch (pełny przebieg)",
                        "Zoom: PWM → sinus po filtrze notch"
                    ))

# Pełny przebieg (próbkowanie współczynnikiem 20 dla czytelności)
dec = 20
fig.add_trace(go.Scatter(x=t_pwm[::dec] * 1000, y=V_A_pwm[::dec],
                        mode='lines', name='PWM (A)',
                        line=dict(color='gray', width=0.6), opacity=0.6),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_pwm[::dec] * 1000, y=V_A_notch[::dec],
                        mode='lines', name='Po filtrze notch (A)',
                        line=dict(color='red', width=2)),
              row=1, col=1)
fig.add_trace(go.Scatter(x=t_pwm[::dec] * 1000, y=V_A_target[::dec],
                        mode='lines', name='Docelowy sinus (A)',
                        line=dict(color='black', width=2, dash='dash')),
              row=1, col=1)

# Zoom
fig.add_trace(go.Scatter(x=t_pwm[zoom_mask] * 1000, y=V_A_pwm[zoom_mask],
                        mode='lines', name='PWM (zoom)',
                        line=dict(color='gray', width=0.6), opacity=0.6),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_pwm[zoom_mask] * 1000, y=V_A_notch[zoom_mask],
                        mode='lines', name='Po filtrze notch (zoom)',
                        line=dict(color='red', width=2)),
              row=2, col=1)
fig.add_trace(go.Scatter(x=t_pwm[zoom_mask] * 1000, y=V_A_target[zoom_mask],
                        mode='lines', name='Docelowy sinus (zoom)',
                        line=dict(color='black', width=2, dash='dash')),
              row=2, col=1)

fig.update_xaxes(title_text='Czas [ms]', row=1, col=1)
fig.update_xaxes(title_text='Czas [ms]', row=2, col=1)
fig.update_yaxes(title_text='Napięcie [V]', row=1, col=1)
fig.update_yaxes(title_text='Napięcie [V]', row=2, col=1)
fig.update_layout(template='plotly_white', height=720, hovermode='x unified')

fig.show()

print(f"Notch: stop_bw = {stop_bw} Hz, harmoniczne 1..{max_harm}")
print(f"RMS po notch:    {rms_notch:.2f} V")
print(f"RMS celu:        {rms_target:.2f} V")
print(f"Maks. odchyłka:  {max_err:.2f} V")

Notch: stop_bw = 2000 Hz, harmoniczne 1..6
RMS po notch:    16.47 V
RMS celu:        15.27 V
Maks. odchyłka:  27.27 V
In [19]:
# Widmo sygnału PWM po filtrze notch vs przed filtrem
import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots

# Używamy tej samej części ustalonej co wcześniej
sig_pwm = V_A_pwm[start_idx:]
sig_notch = V_A_notch[start_idx:]
Nsig = len(sig_pwm)
freqs = np.fft.fftfreq(Nsig, dt_pwm)
pos = freqs > 0
freq_pos = freqs[pos]
fft_pwm = np.abs(np.fft.fft(sig_pwm)[pos]) * 2 / Nsig
fft_notch = np.abs(np.fft.fft(sig_notch)[pos]) * 2 / Nsig

# Poziom głównej harmonicznej PWM przed/po
idx_pwm = np.argmin(np.abs(freq_pos - f_pwm))
amp_pwm = fft_pwm[idx_pwm]
amp_notch = fft_notch[idx_pwm]
atten_db = 20 * np.log10(amp_notch / amp_pwm) if amp_pwm > 0 else np.nan

fig = make_subplots(rows=2, cols=1, vertical_spacing=0.12,
                    subplot_titles=("Widmo: przed vs po filtrze notch (pełny zakres)",
                                    "Zoom 0-5 kHz"))

# Pełny zakres (log y)
fig.add_trace(go.Scatter(x=freq_pos/1000, y=fft_pwm,
                         mode='lines', name='PWM przed', line=dict(color='gray', width=1)),
              row=1, col=1)
fig.add_trace(go.Scatter(x=freq_pos/1000, y=fft_notch,
                         mode='lines', name='PWM po notch', line=dict(color='red', width=1.3)),
              row=1, col=1)
fig.update_yaxes(type='log', row=1, col=1)
fig.update_xaxes(title_text='Częstotliwość [kHz]', row=1, col=1)
fig.update_yaxes(title_text='Amplituda [V]', row=1, col=1)
fig.add_vline(x=f_electrical/1000, line_dash='dash', line_color='green', row=1, col=1)
fig.add_vline(x=f_pwm/1000, line_dash='dash', line_color='blue', row=1, col=1)

# Zoom niskie częstotliwości
zoom_mask = freq_pos <= 5000
fig.add_trace(go.Scatter(x=freq_pos[zoom_mask], y=fft_pwm[zoom_mask],
                         mode='lines', name='PWM przed (zoom)', line=dict(color='gray', width=1)),
              row=2, col=1)
fig.add_trace(go.Scatter(x=freq_pos[zoom_mask], y=fft_notch[zoom_mask],
                         mode='lines', name='PWM po notch (zoom)', line=dict(color='red', width=1.5)),
              row=2, col=1)
fig.update_xaxes(title_text='Częstotliwość [Hz]', row=2, col=1)
fig.update_yaxes(title_text='Amplituda [V]', row=2, col=1)

fig.update_layout(template='plotly_white', height=720, hovermode='x unified')
fig.show()

print(f"Harmoniczna f_PWM: przed = {amp_pwm:.3f} V, po = {amp_notch:.3f} V, tłumienie = {atten_db:.1f} dB")

Harmoniczna f_PWM: przed = 17.075 V, po = 0.000 V, tłumienie = -101.4 dB

próba usuniecia harmonicznych filtrem LP powyżej 30kHz¶

In [20]:
# Filtr FFT: wycięcie powyżej 30 kHz (prostokątny low-pass) + wizualizacja czasu i widma
import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots

fc_lp_fft = 30000  # Hz, granica odcięcia

# FFT na pełnym sygnale (aby mieć V_A_lp30 tej samej długości co V_A_pwm)
N_full = len(V_A_pwm)
freqs_full = np.fft.fftfreq(N_full, dt_pwm)
fft_vals_full = np.fft.fft(V_A_pwm)

# Maska LP: zostawiamy |f| <= fc_lp_fft
mask_lp = np.abs(freqs_full) <= fc_lp_fft
fft_lp_full = fft_vals_full * mask_lp
V_A_lp30 = np.fft.ifft(fft_lp_full).real

# Metryki czasu
rms_lp = np.sqrt(np.mean(V_A_lp30**2))
rms_target = np.sqrt(np.mean(V_A_target**2))
max_err = np.max(np.abs(V_A_lp30 - V_A_target))

# === WIZUALIZACJA W DZIEDZINIE CZASU ===
zoom_start = 0.005
zoom_end = 0.008
zoom_mask_time = (t_pwm >= zoom_start) & (t_pwm <= zoom_end)

fig_time = make_subplots(rows=2, cols=1, vertical_spacing=0.12,
                         subplot_titles=(
                             "PWM vs filtr LP 30 kHz (pełny przebieg)",
                             "Zoom: PWM → sinus po filtrze LP 30k"
                         ))

# Pełny przebieg (próbkowanie współczynnikiem 20 dla czytelności)
dec = 20
fig_time.add_trace(go.Scatter(x=t_pwm[::dec] * 1000, y=V_A_pwm[::dec],
                              mode='lines', name='PWM (A)',
                              line=dict(color='gray', width=0.6), opacity=0.6),
                   row=1, col=1)
fig_time.add_trace(go.Scatter(x=t_pwm[::dec] * 1000, y=V_A_lp30[::dec],
                              mode='lines', name='Po filtrze LP 30k (A)',
                              line=dict(color='red', width=2)),
                   row=1, col=1)
fig_time.add_trace(go.Scatter(x=t_pwm[::dec] * 1000, y=V_A_target[::dec],
                              mode='lines', name='Docelowy sinus (A)',
                              line=dict(color='black', width=2, dash='dash')),
                   row=1, col=1)

# Zoom
fig_time.add_trace(go.Scatter(x=t_pwm[zoom_mask_time] * 1000, y=V_A_pwm[zoom_mask_time],
                              mode='lines', name='PWM (zoom)',
                              line=dict(color='gray', width=0.6), opacity=0.6),
                   row=2, col=1)
fig_time.add_trace(go.Scatter(x=t_pwm[zoom_mask_time] * 1000, y=V_A_lp30[zoom_mask_time],
                              mode='lines', name='Po filtrze LP 30k (zoom)',
                              line=dict(color='red', width=2)),
                   row=2, col=1)
fig_time.add_trace(go.Scatter(x=t_pwm[zoom_mask_time] * 1000, y=V_A_target[zoom_mask_time],
                              mode='lines', name='Docelowy sinus (zoom)',
                              line=dict(color='black', width=2, dash='dash')),
                   row=2, col=1)

fig_time.update_xaxes(title_text='Czas [ms]', row=1, col=1)
fig_time.update_xaxes(title_text='Czas [ms]', row=2, col=1)
fig_time.update_yaxes(title_text='Napięcie [V]', row=1, col=1)
fig_time.update_yaxes(title_text='Napięcie [V]', row=2, col=1)
fig_time.update_layout(template='plotly_white', height=720, hovermode='x unified')

fig_time.show()

print(f"LP FFT 30 kHz: stop > {fc_lp_fft/1000} kHz")
print(f"RMS po LP:       {rms_lp:.2f} V")
print(f"RMS celu:        {rms_target:.2f} V")
print(f"Maks. odchyłka:  {max_err:.2f} V")

# === ANALIZA WIDMOWA ===
# Część ustalona sygnału dla FFT
sig_pwm = V_A_pwm[start_idx:]
sig_lp = V_A_lp30[start_idx:]
Nsig = len(sig_pwm)
freqs = np.fft.fftfreq(Nsig, dt_pwm)
fft_vals = np.fft.fft(sig_pwm)
fft_lp = np.fft.fft(sig_lp)

# Widma (tylko dodatnie częstotliwości)
pos = freqs > 0
freq_pos = freqs[pos]
fft_pwm_abs = np.abs(fft_vals[pos]) * 2 / Nsig
fft_lp_abs = np.abs(fft_lp[pos]) * 2 / Nsig

# TÅ‚umienie na f_pwm
idx_pwm = np.argmin(np.abs(freq_pos - f_pwm))
amp_pwm = fft_pwm_abs[idx_pwm]
amp_lp = fft_lp_abs[idx_pwm]
atten_db = 20 * np.log10(amp_lp / amp_pwm) if amp_pwm > 0 else np.nan

fig_freq = make_subplots(rows=2, cols=1, vertical_spacing=0.12,
                         subplot_titles=("Widmo: przed vs LP 30 kHz (pełny zakres)",
                                        "Zoom 0-5 kHz"))

# Pełny zakres (log y)
fig_freq.add_trace(go.Scatter(x=freq_pos/1000, y=fft_pwm_abs,
                              mode='lines', name='PWM przed', line=dict(color='gray', width=1)),
                   row=1, col=1)
fig_freq.add_trace(go.Scatter(x=freq_pos/1000, y=fft_lp_abs,
                              mode='lines', name='PWM po LP 30k', line=dict(color='red', width=1.3)),
                   row=1, col=1)
fig_freq.update_yaxes(type='log', row=1, col=1)
fig_freq.update_xaxes(title_text='Częstotliwość [kHz]', row=1, col=1)
fig_freq.update_yaxes(title_text='Amplituda [V]', row=1, col=1)
fig_freq.add_vline(x=f_electrical/1000, line_dash='dash', line_color='green', row=1, col=1)
fig_freq.add_vline(x=f_pwm/1000, line_dash='dash', line_color='blue', row=1, col=1)
fig_freq.add_vline(x=fc_lp_fft/1000, line_dash='dot', line_color='black', row=1, col=1,
                   annotation_text='fc=30 kHz')

# Zoom niskie częstotliwości
zoom_mask_freq = freq_pos <= 5000
fig_freq.add_trace(go.Scatter(x=freq_pos[zoom_mask_freq], y=fft_pwm_abs[zoom_mask_freq],
                              mode='lines', name='PWM przed (zoom)', line=dict(color='gray', width=1)),
                   row=2, col=1)
fig_freq.add_trace(go.Scatter(x=freq_pos[zoom_mask_freq], y=fft_lp_abs[zoom_mask_freq],
                              mode='lines', name='PWM po LP 30k (zoom)', line=dict(color='red', width=1.5)),
                   row=2, col=1)
fig_freq.update_xaxes(title_text='Częstotliwość [Hz]', row=2, col=1)
fig_freq.update_yaxes(title_text='Amplituda [V]', row=2, col=1)

fig_freq.update_layout(template='plotly_white', height=720, hovermode='x unified')
fig_freq.show()

print(f"\nWIDMO - Harmoniczna f_PWM: przed = {amp_pwm:.3f} V, po = {amp_lp:.3f} V, tłumienie = {atten_db:.1f} dB")

LP FFT 30 kHz: stop > 30.0 kHz
RMS po LP:       20.53 V
RMS celu:        15.27 V
Maks. odchyłka:  30.91 V

WIDMO - Harmoniczna f_PWM: przed = 17.075 V, po = 17.075 V, tłumienie = -0.0 dB