Analysis of Sinusoidal Pulse Width Modulation of AC Signal

Hello friends! Today, we will analyze sinusoidal pulse-width modulation of an AC signal in MATLAB. We will generate a carrier, compare it with a reference, calculate a switched voltage and examine the resulting current and frequency spectrum. The aim is to understand what the equations mean rather than simply obtain a collection of plots.

The original tutorial gates an existing sinusoidal supply on and off. That is an idealized AC-chopper model, which differs from a DC-fed inverter that synthesizes AC by switching a DC link. We will preserve that original idea and make the switching assumptions explicit so the calculated voltage, current and power factor describe one consistent circuit model.

We will also correct several common analysis mistakes: confusing pulse count with period, treating per-unit quantities as limited to one, assuming the first FFT bin is always the fundamental, and calling THD a threshold. A clearly defined waveform and measurement window are essential before any reported result can be trusted.

Define the Switching Model First

Let the input be v_s(t) = √2V_s,rmssin(2πf₀t). A binary gate g(t) selects either the source voltage or a zero-voltage freewheeling state at the load. Our model therefore uses v_o(t) = g(t)v_s(t).

The load is a series resistor and inductor. During the zero-voltage state, current continues circulating in an ideal freewheeling path. The source current is i_s(t) = g(t)i_L(t). This requires an ideal bidirectional switching arrangement; it is not a claim that one transistor and one diode can perform all the necessary AC commutation.

We neglect switch losses, dead time, diode drops, parasitic elements and source impedance. Those simplifications are acceptable for studying the algorithm if they remain visible. A physical AC chopper needs a separately designed power stage and commutation strategy.

Generate a Carrier and Reference

We use a triangular carrier between zero and one and compare it with m_a|sin(2πf₀t)|. The gate is high when the reference is at least as large as the carrier. Increasing the reference widens the pulses, while the carrier frequency determines how often switching opportunities occur.

For this definition, keep m_a between zero and one. That is a constraint on this modulation scheme, not a universal restriction on every quantity called a modulation index. Also distinguish the carrier from a sawtooth: a triangle has rising and falling ramps, whereas a sawtooth returns abruptly after one ramp.

If f₀ is 50 Hz and f_c is 2,000 Hz, there are 40 carrier cycles per fundamental period and 20 per half-period. The fundamental period remains 0.02 seconds. A pulse count of 40 is dimensionless and must not be described as a 40-second period.

The original input dialog is shown below. It requests frequency, modulation index, load phase angle and pulse count. The revised code uses explicit parameter assignments so a run can be reproduced without remembering a series of interactive answers.

A Sinusoidal Reference Does Not Guarantee a Sinusoidal Output

In our chopper, the local average duty ratio is approximately D(t) = m_a|sin(ωt)| when switching is sufficiently fast. Multiplying that duty by the already sinusoidal source gives an approximate low-frequency waveform proportional to sin(ωt)|sin(ωt)|. That shape is not a sine wave.

This matters because the familiar DC-link inverter relation between modulation index and fundamental amplitude cannot simply be copied into this model. Under the high-carrier-frequency averaging approximation, the fundamental peak of this chopper waveform is 8m_aV_s,peak/(3π). Other low-order components remain.

A constant duty ratio applied to the same ideal AC source would instead give an average waveform proportional to that source sine. Whether that is the appropriate strategy depends on the application. The lesson is to derive the switched output from the actual topology and switching rule before judging its spectrum.

Define the Load with Consistent Units

For a series RL load at the fundamental frequency, Z = R + jω₀L. If you specify impedance magnitude |Z| and phase angle φ, then R = |Z|cos(φ) and L = |Z|sin(φ)/(2πf₀). The phase must be converted to radians when using MATLAB's sin and cos functions.

For |Z| = 10 ohms, φ = 60 degrees and f₀ = 50 Hz, R = 5 ohms and L is approximately 0.02757 H. The time constant is L/R ≈ 5.51 milliseconds. The impedance angle describes the load at 50 Hz; harmonics encounter different inductive reactance.

MATLAB's pi is approximately 3.14159. It is not 2.14. Likewise, per-unit notation is a ratio to a selected base, not a number that must remain below one. A voltage of 1.1 per unit simply means 110% of the chosen voltage base.

This example uses volts, amperes, ohms and henries directly. If using per-unit values, define voltage, current, impedance and frequency bases first. A numerical inductance derived from an unexplained per-unit impedance cannot automatically be labeled in henries.

Calculate Current Through the RL Load

The load equation is L di/dt + Ri = v_o. Assuming the applied voltage remains constant over one numerical sample interval Δt, the exact first-order update over that interval is:

i[k+1] = ai[k] + (1 − a)v_o[k]/R, with a = e^−RΔt/L.

The update is exact for the assumed constant voltage during each interval, not for the continuously changing source between samples. Reducing Δt improves the representation of both the source and switching instants. It also reduces pulse-width quantization error.

Start with a defined current, then allow transients to decay before measuring a steady-state spectrum. Choosing an arbitrary negative initial current and analyzing immediately can make the result mostly describe startup energy rather than periodic operation. Our code uses enough cycles to exceed twelve load time constants before the final analysis window.

A purely resistive load requires its own branch, i = v/R, because the exponential formula divides by L. The example below requires R and L to be positive; it does not silently accept an invalid phase angle or zero inductance.

Complete MATLAB Analysis Code

Save the following as acChopperAnalysis.m and call results = acChopperAnalysis. It uses a coherent integer number of carrier cycles per fundamental cycle and analyzes the final complete fundamental period. No Signal Processing Toolbox function is required by this example.

function results = acChopperAnalysis
    f0 = 50;
    carrierRatio = 40;
    samplesPerCarrier = 200;
    ma = 0.8;
    sourceRms = 24;
    impedanceMagnitude = 10;
    phi = 60*pi/180;

    R = impedanceMagnitude*cos(phi);
    L = impedanceMagnitude*sin(phi)/(2*pi*f0);
    assert(f0 > 0 && sourceRms > 0 && R > 0 && L > 0);
    assert(ma > 0 && ma <= 1);
    assert(carrierRatio == floor(carrierRatio));
    assert(samplesPerCarrier == floor(samplesPerCarrier));
    assert(carrierRatio > 0 && samplesPerCarrier >= 20);

    N = carrierRatio*samplesPerCarrier;
    fs = f0*N;
    dt = 1/fs;
    tau = L/R;
    cycles = max(20,ceil(12*tau*f0)+2);
    count = cycles*N;
    k = (0:count-1)';
    t = k*dt;
    theta = 2*pi*mod(k,N)/N;
    carrierPhase = mod(k,samplesPerCarrier)/samplesPerCarrier;
    carrier = 1-abs(2*carrierPhase-1);
    reference = ma*abs(sin(theta));
    gate = reference >= carrier;
    vs = sqrt(2)*sourceRms*sin(theta);
    vo = double(gate).*vs;

    iLoad = zeros(count,1);
    a = exp(-R*dt/L);
    b = -expm1(-R*dt/L)/R;
    for n = 2:count
        iLoad(n) = a*iLoad(n-1)+b*vo(n-1);
    end
    iSource = double(gate).*iLoad;
    take = count-N+1:count;
    v = vs(take); out = vo(take);
    il = iLoad(take); source = iSource(take);

    [Vrms,V1,Vdc,Vthd] = spectrumMetrics(out);
    [Irms,I1,Idc,Ithd] = spectrumMetrics(il);
    [Isrms,Is1,Isdc,Isthd] = spectrumMetrics(source);
    inputRms = sqrt(mean(v.^2));
    realPower = mean(v.*source);
    PF = realPower/(inputRms*Isrms);
    Fv = fft(v)/N;
    Fi = fft(source)/N;
    displacement = angle(Fv(2)*conj(Fi(2)));

    results = table([Vrms;Irms;Isrms],[V1;I1;Is1], ...
        [Vdc;Idc;Isdc],100*[Vthd;Ithd;Isthd], ...
        'VariableNames',{'RMS','FundamentalRMS','DC','THDpercent'}, ...
        'RowNames',{'OutputVoltage','LoadCurrent','SourceCurrent'});
    disp(results);
    fprintf('Input power: %.6g W\nPower factor: %.6g\n',realPower,PF);
    fprintf('Fundamental displacement: %.6g degrees\n', ...
        displacement*180/pi);
    localTime = (0:N-1)'*dt;
    subplot(3,1,1); plot(localTime,carrier(take), ...
        localTime,reference(take)); ylabel('Carrier / reference');
    subplot(3,1,2); plot(localTime,v,localTime,out); ylabel('Voltage (V)');
    subplot(3,1,3); plot(localTime,il,localTime,source);
    ylabel('Current (A)'); xlabel('Time (s)');
end

function [total,fundamental,dc,distortion] = spectrumMetrics(x)
    x = x(:);
    N = numel(x);
    F = fft(x)/N;
    total = sqrt(mean(x.^2));
    dc = real(F(1));
    fundamental = sqrt(2)*abs(F(2));
    if fundamental <= 1e-12
        distortion = NaN;
    else
        remainder = max(0,total^2-dc^2-fundamental^2);
        distortion = sqrt(remainder)/fundamental;
    end
end

The code retains the input, output and source-current signals separately. It uses the gate state directly to identify source conduction rather than testing whether output voltage happens to equal zero. At an AC zero crossing, voltage can be zero even while a switch is commanded on.

The absolute tolerance in the helper prevents division by an effectively zero fundamental for this small numerical example. A reusable measurement library should choose tolerances relative to the signal scale and also validate every parameter as finite. The source above is an educational calculation, not gate-drive firmware.

Read RMS and FFT Results Correctly

For uniformly sampled data, total RMS is √mean(x²). It includes DC, the fundamental and all represented harmonics. The fundamental RMS describes only the sinusoidal component at f₀. A heavily switched waveform can have a large total RMS without having a correspondingly large useful fundamental.

MATLAB's first FFT array element is DC. In this code, the second element corresponds to 50 Hz only because the analysis window contains exactly one 50 Hz period. If the window contained ten periods, the fundamental would be at index 11. Never discard DC and assume the next element is the fundamental for every record length.

For a real coherent signal, the positive-frequency fundamental has a matching negative-frequency component. Its peak amplitude is 2|F(2)| after division by sample count, and RMS is √2|F(2)|. The DC and Nyquist terms require different treatment. See the official FFT amplitude-scaling documentation.

THD means total harmonic distortion. In a coherent periodic record, THD = √(Σ_h≥2X_h,rms²)/X_1,rms. Our helper obtains the remainder from total mean-square after subtracting DC and fundamental contributions. This includes all remaining energy represented in the sampled record.

With noise, nonperiodic transients or interharmonics, that remainder is not strictly harmonic-only distortion. State the assumptions and bandwidth instead of calling every residual-energy ratio pure THD. Also distinguish a ratio, percentage and decibels: a ratio of 0.1 is 10%, or −20 dB relative to the fundamental amplitude.

MATLAB's thd function has its own harmonic-selection, windowing and output-unit conventions. It should not be expected to reproduce an all-bin calculation unless the measurement definitions are aligned.

Calculate True Input Power Factor

True power factor is PF = mean(v_si_s)/(V_s,rmsI_s,rms). This compares real input power with apparent power. Calculate it from matching samples of voltage and source current over the same steady-state interval.

For a sinusoidal supply and suitable periodic current, PF can also be expressed as displacement factor multiplied by current distortion factor: cos(φ₁)I_1,rms/I_rms. The load impedance angle alone does not determine the input power factor of a switched converter.

Use the relative complex phase of the fundamental voltage and current. A one-argument arctangent can lose quadrant information. The code uses the angle of F_vconj(F_i), which also makes the phase reference explicit.

Review the Original Plots and Improve the Experiment

The historical waveform plot shows the original source, reference and switching construction. It remains useful for visual comparison, but it was produced by the earlier algorithm rather than the revised code above.

The next figure combines voltage, current and spectral information. Read the units and harmonic axes carefully: a plotted FFT array index is not automatically a frequency in hertz or a harmonic order unless the conversion is defined.

The original numerical output is retained below. Treat it as a historical run with its own assumptions, not as independently verified results for the revised parameters.

Repeat the calculation with 200 and 400 samples per carrier cycle while keeping the physical parameters fixed. Compare RMS, fundamental amplitude, THD and power factor. If the results change materially, the coarser simulation has not resolved the pulses adequately.

Check energy consistency as well. In periodic steady state, an ideal inductor has zero net stored-energy change over a cycle, so input real power should approach R mean(i_L²) as numerical resolution improves. A significant unexplained mismatch points to inconsistent switching paths, transient contamination or sampling error.

Choose a Meaningful Frequency Range

The sample frequency determines the Nyquist limit, while the record duration determines FFT bin spacing. In our example, increasing samples per carrier raises the highest represented frequency, but a one-period record still has bins separated by 50 Hz. These are different improvements: more samples resolve switching edges, whereas a longer observation interval provides finer frequency spacing.

Zero-padding makes a plotted spectrum appear more densely sampled between existing frequency points, but it does not create additional observed cycles or separate components that the original record cannot resolve. If you extend the record, update the fundamental-bin selection instead of retaining the hardcoded second FFT element.

A coherent record contains an integer number of periods of the repeating waveform. Here, the integer carrier ratio and repeating discrete construction help maintain that condition. With a noninteger carrier ratio, the full switching pattern may repeat over several source cycles or may not fit the selected window at all. In that case, a simple single-cycle harmonic interpretation becomes inappropriate.

Measure Pulse Width Without Guessing Crossing Points

The older code searched for approximate equality between carrier and reference samples. That approach can miss a crossing when no sample lies within the chosen tolerance, or record several neighboring samples as separate matches. The binary gate sequence provides a clearer starting point: find its rising and falling transitions and pair them in time order.

A discrete pulse width is the number of high sample intervals multiplied by Δt. Pulses touching the beginning or end of a record need special handling because their other edge may lie outside the window. Do not silently pair the first falling edge with a later rising edge and produce a negative duration.

For a 2,000 Hz carrier with 200 samples per cycle, Δt is 2.5 microseconds. An edge located only to the nearest sample carries that time quantization into the pulse width. Interpolating the continuous carrier-reference crossing can improve an event estimate, but only if the interpolated waveform assumptions are stated and the switching simulation uses them consistently.

Check Limiting Cases Before Trusting a Sweep

With a manually forced gate of one, the model should reduce to a sinusoidal source feeding the RL load. Its steady-state current RMS should approach V_s,rms/|Z| and its phase lag should approach the selected impedance angle. This is a useful independent check of the current update and phase calculation.

With a manually forced gate of zero, the source current and applied load voltage should be zero, while any initial load current decays through the modeled freewheeling loop. These forced-gate cases are diagnostic modifications, not the ordinary ma settings of the displayed code.

Finally, compare two runs using different initial currents after adequate settling time. Their final periodic measurements should converge. If they do not, investigate the settling interval, numerical update or switching construction before presenting a parameter sweep as a physical conclusion.

Useful Extensions and Related Tutorials

Try changing modulation index, carrier ratio and load phase separately. Record each parameter set alongside its metrics. Then compare this gated-AC model with a separately defined DC-link inverter, where the available switch voltages and modulation equations differ.

You can build a parameter interface with our MATLAB GUI tutorial or explore serial communication. Other signal-processing examples include correlation-based speech recognition and DTMF decoding.

The articles on hexapod simulation, P-code and Simulink packaging provide further MATLAB project context. Keep numerical assumptions documented when sharing any analysis.

When reporting a spectrum, include the source frequency, carrier frequency, sample frequency, analysis duration and harmonic range. Also state whether the waveform is measured, simulated with switching devices or generated algebraically. These details allow another reader to reproduce the calculation and explain differences without assuming that one graph is wrong simply because it looks less smooth.

Frequently Asked Questions

Why is the output not a pure sine wave?

It contains switching pulses, and this particular sinusoidal-duty AC-chopper rule also creates a nonsinusoidal low-frequency average. Filtering cannot correct an incorrectly assumed modulation law.

Can per-unit voltage exceed one?

Yes. Per-unit means a ratio to a defined base. Values above one represent quantities above that base.

Why must startup samples be excluded?

The initial inductor current introduces a transient that is not part of the periodic steady-state waveform. Including it changes RMS and spectral measurements.

Does this code define a safe power circuit?

No. It defines an ideal numerical switching model. Device selection, bidirectional commutation, protection and gate timing require a separate hardware design.

Were the revised results executed in MATLAB?

No fresh MATLAB run is claimed here. The source and measurement definitions are provided for reproduction; the retained screenshots belong to the earlier tutorial.