PhD Course on  Harmonics in Power Electronics and Power Systems

This course provides a broad overview of power system harmonic problems, methods of analyzing, measuring and effectively mitigating them. Several extended simulation and data processing tools, among others DIgSILENT PowerFactory, Matlab/Simulink or LabVIEW are used to assess and study the harmonic distortion at different points of power networks.
The results of analytical investigation and simulations are validated against measurements applying sophisticated data processing techniques. Furthermore, deep understanding of hardware considerations regarding har- monic measurements in harsh industrial environment is given, using specialized equipment, for in- stance GPS-synchronized measuring instruments.

The course covers the following topics:

  • Power Quality definitions. Generation mechanism of power system harmonics. Harmonic indices.
  • Voltage vs. current distortion as well as parallel vs. series resonance in modern power systems. Point of Common Coupling (PCC).
  • Sources and effects of harmonic distortion.
  • Harmonic measuring instruments and measuring procedures in LV, MV and HV networks.
  • Mathematical tools and theories for analyzing distorted waveforms. Signal processing and uncertainty analysis.
  • Modelling of classical power system components. Harmonic analysis.
  • Modelling of grid-connected converters for harmonic analysis purposes and their application in modern power systems including e.g. offshore wind power plants.
  • Harmonic load-flow, frequency scan and time domain simulations. Linear and nonlinear analysis techniques.
  • Steady-state harmonics vs. harmonic stability. Small-signal representation, sequence and frequency coupling.
  • Software tools for harmonic analysis.
  • Precautionary (preventive) and corrective (remedial) harmonic mitigation techniques. Passive and active line filters. Filter design.

Organizer: Professor Claus Leth Bak
Lecturers: Christian Frank Flytkjær from Energinet and Łukasz H. Kocewiak from Ørsted

Harmonic current of 6-pulse rectifier supplying a resistive load
Figure 1 Harmonic current of 6-pulse rectifier supplying a resistive load.

Oversampling (resampling) data in harmonic data processing

When using a rectangular window, it is important to synchronize the measurement window accurately with the power system frequency to achieve integer multiple of periods in analysed time series (i.e. window period is and integer multiple of analysed frequency component natural period). For example, if the power system frequency is 50.2 Hz whereas the window size is 200 ms, the fundamental frequency spectral line of the discrete Fourier transform is no longer projected (represented) by one complex vector in the orthogonal basis but spanned over the whole basis. One can even say that the power system frequency in the estimated spectrum becomes an interharmonic with spectral leakage as a consequence. Therefore appropriate data should be prepared before the projection into frequency domain.

In order to limit possible data processing errors the signal should be adjusted before spectral analysis. Power system frequency variation can be limited by application of various interpolation methods. It is strongly recommended to always oversample the analysed signal. In this particular case acquired waveforms with sample rate of 44.1 kS/s/ch are oversampled up to 51.2 kS/s/ch. According to this procedure 1024 samples per cycle is obtained, which is an integer power of 2 and can be used to apply fast Fourier transform. Please note that for 10-cycle data blocks discrete Fourier transform should be applied.

Frequency estimation
Figure 1 Single tone frequency. (a) Measured power system frequency variation. (b) Frequency variation after resampling (spline interpolation).

Oversampling of each 10-cycle data blocks before spectrum estimation significantly limits fundamental frequency variation over analysed window and therefore the stationarity assumption is strengthened for frequency component of interest (i.e. linked with the fundamental frequency). As it can be seen in Figure 1 resampling with high precision definitely improves signal quality. The figure presents analysed one minute of acquired voltage waveform and each point is an estimated fundamental frequency over 10-cycle rectangular window. The discrepancy as a difference between two measured values, given in Figure 1(b), is acceptably small and the sample standard deviation is 4.3548•10-6. Please note that scale 0.2s used on abscissa (i.e. horizontal axis) in Figure 1 denotes the step size between neighbouring samples. If there are 300 samples in presented results (as in Figure 1) the total duration is 1 minute (i.e. 0.2s∗300=60s).

It can be seen that both factors such as main tone detection algorithm as well as interpolation algorithm are crucial in appropriate harmonic magnitude and phase estimation. Various interpolation techniques can give different results. The most sophisticated unfortunately are usually the most time consuming. Therefore taking into consideration enormous amount of data to process from few months of multipoint measurements also the interpolation should be optimized and agreement between acceptable accuracy and calculation speed should be found.

Roughly speaking, during interpolation process new samples between existing samples are computed. Different methods were used in data processing such as nearest (coerce) method, linear method, spline method, cubic method, and finite impulse response filter method [1]. The nearest method finds the point nearest to xi  in an X array and then assigns the corresponding y value in an interpolated Y to yi. The linear method interpolates yi on the line which connects the two points (xj,xj+1) when the corresponding location xi of interpolated value yi which is located between the two points in an X array. The spline method known also as cubic spline method is almost always the most suitable solution. Within the method the third-order polynomial for each interval between two adjacent points in X (xj,xj+1) is derived. The polynomials have to meet the following conditions: the first and second derivatives at xj and xj+1 are continuous, the polynomials pass all the specified data points, the second derivatives at the beginning and end are zero. The cubic Hermitian spline method is the piecewise cubic Hermitian interpolation. This method is similar to cubic spline interpolation and derives a third-order polynomial but in Hermitian form for each interval (xj,xj+1) and ensures only the first derivatives of the interpolation polynomials are continuous. Compared to the cubic spline method, the cubic Hermitian method is characterized by better local property. The cubic Hermite interpolation method guarantees that the first derivative of the interpolant is continuous and sets the derivative at the endpoints in order to preserve the original shape and monotonicity of the data [2]. Also interpolation based on finite impulse response filter is popular [3], [4]. This method provides excellent results also for frequency analysis, although it is more intensive computationally [2].

Order Interpolation method Computation time [ms] Marker
1 Linear 110 x
2 Coerce 70 *
3 Cubic spline 130
4 Hermitian spline 430 +
5 FIR filter 290

Table above shows the computation time of data interpolation with different algorithms. In each of the cases the presented time is measured for resampling of 10 cycles of measured data. It can be seen that the simplest algorithms exhibit less computation burden. On the other hand in many cases simple algorithms cannot give satisfactory interpolation results and thus affect the harmonic calculation results. Figure 2 shows how different interpolation techniques affect harmonic magnitude estimation from measurements carried out at the LV side of the wind turbine transformer. It can be seen that the coerce method is affected the worst while other algorithms give comparable results.

Oversampled measurements
Figure 2 Different interpolation techniques used in oversampling of waveforms measured at LV side of the wind turbine transformer.

It was observed that in case of more distorted waveforms (e.g. gird-side converter output voltage) interpolation method choice plays a crucial role. In Figure 3 one can see that also linear interpolation and cubic Hermitian spline interpolation do not give satisfactory results. Therefore the most suitable for harmonic components estimation appear to be cubic spline interpolation as well as finite impulse response filter interpolation. Since the spline method (as presented in table above) is less time consuming, this method was used in further data processing. The selected measurement period was selected when the power system frequency were varying below 50 Hz or above 50 Hz.

Oversampled measurements
Figure 3 Different interpolation techniques used in oversampling of waveforms measured at the grid-side converter AC terminals.

[1] S. C. Chapra and R. Canale, Numerical Methods for Engineers: With Software and Programming Applications, 6th ed. McGraw-Hill Science, 2009.
[2] National Instruments, "LabVIEW 2011 Help," National Instruments Manual, 2011.
[3] R. A. Losada, "Digital Filters with MATLAB," The Mathworks, 2008.
[4] R. W. Schafer and L. R. Rabiner, "A digital signal processing approach to interpolation," Proceedings of the IEEE, vol. 61, no. 6, pp. 692-702, Jun. 1973.


Commercial power quality meters

The use of a rectangular window requires that the measurement window is synchronized with the actual power system frequency, hence the use of a 10-cycle window instead of a window of exactly  200 ms. The IEC standard [1] requires that 10 cycles correspond with an integer number of samples within 0.03%. To ensure synchronism between the measurement window and the power system frequency, most power quality meters use a phase-locked loop generating a sampling frequency that is an integer multiple of the actual power system frequency. One of the commonly used in the company power quality meters is Elspec G4500 which provides full functionality regarding measurements of power quality. An exemplary connection of such measurement equipment in there phase system is presented in Figure 1 [2].

Elspec G4500 connection
Figure 1 Exemplary connection of Elspec equipment in wind turbines.

Before power quality indices are calculated by the Elspec equipment, acquired data is processed. Processing stage comprises of fundamental frequency detection, sample rate adjustment according to the detected frequency, and Fourier decomposition applied. The approach of sample rate adjustment is to adjust the finite orthogonal Hilbert basis in order to express each of frequency components in the Fourier space only by one vector from the Hilbert basis. The whole data acquisition, processing and logging process is briefly presented in Figure 2.

A synchronization (i.e. sample rate adjustment according to the fundamental frequency) error leads to cross-talk between different harmonic frequencies. The 50 hz component is by far the dominating component in most cases so that the main concern is the cross-talk from the 50 hz component to higher order components. From the other side, resampling process can affect frequency components which are not integer multiple of the fundamental frequency.

This phenomenon can be easily seen in the spectrum of pulse width modulated voltage source converters with fixed frequency ratio. Since generated output voltage of the wind turbine is a function of the fundamental frequency (fo) and the carrier frequency (fc), results obtained by the Elspec system are incorrect to some extent. The magnitudes of Fourier transform harmonic components are the more inaccurate the higher significant is the fundamental variation. Even if the result of applying the discrete Fourier transform to the basic window is a spectrum with 5-hz spacing between frequency components and the spectrum thus contains both harmonics and interharmonics, the results can be sometimes significantly inaccurate.

Power quality meter data processing
Figure 2 Data processing algorithm used by Elspec equipment.

As a good example of this is one of the most significant sideband harmonic components from the first carrier group. Exemplary results of measurements form the LV side of the wind turbine transformer can be seen in Figure 3. As it was mentioned previously the wind turbine frequency ratio is mf=49 and the analyzed sideband harmonic component is of frequency fc+2fo. Three scatter plots present the same harmonic component measured during the same period using different data acquisition devices and processing techniques.

Sideband harmonic - not resampled signal Sideband harmonic - resampled signal
(a) (b)
Sideband harmonic - Elspec G4500

Figure 3 Sideband harmonic component affected by different processing techniques: (a) sideband harmonic calculated in post-processing from resampled signal, (b) sideband harmonic calculated in post-processing from original signal, (c) sideband harmonic calculated on-line by power quality meter.

Data processing results presented in Figure 3 show how easily inappropriately applied processing techniques can provide wrong results. Results from Figure 3(a) show sideband harmonic component measured using results obtained from the measurement campaign with direct Fourier decomposition (i.e. without sample rate adjustment). From Figure 3(b) presents the same data but resampled and later discrete Fourier transform is applied, Figure 3(c) describe sideband harmonic component magnitude obtained by the Elspec measurement system.

It can be observed that results using various processing approach provide different results. Due to the fact that the frequency of sideband harmonic components generated by modulation with constant carrier frequency does not vary significantly and Fourier decomposition without earlier sample rate adjustment gives the most appropriate results. This is presented in Figure 3(a) and only small magnitude variation affected by nonlinear relation between the modulation index () and the sideband harmonic components as well as measurement and data processing (i.e. small spectral leakage) errors. Completely different and unacceptable results are seen in Figure 3(b) where analysed waveform is resampled. One can observe that due to significant spectral leakage the estimated magnitude sometimes can be even equal to zero. Therefore sometimes power quality meters can provide values significantly affected by processing errors. This behaviour is present in the scatter plot from Elspec measurements (Figure 3(c)). It is important to emphasize that the algorithm in the power quality meter applies lossy compression which also determines estimated magnitudes. Estimated harmonic components are assumed to be insignificant and not saved (i.e. set to zero) if the magnitude is lower than a certain threshold which is defined based on measured waveform distortion and maximum allowed database storage capacity per month. Such limitations provides scatter plots as in Figure 3(c) which is similar to Figure 3(b) but modified due to averaging and magnitudes below the threshold artificially set to zero.

[1] "Electromagnetic compatibility (EMC) - Part 4-30: Testing and measurement techniques - Power quality measurement methods," International Electrotechnical Commission Standard IEC 61000-4-30, 2008.
[2] P. Nisenblat, A. M. Broshi, and O. Efrati, "Methods of compressing values of a monitored electrical power signal," U.S. Patent 7,415,370 B2, Aug. 19, 2008.


Shapes in deep shadows and high lights in Matlab

Normally I work in different toolboxes of Matlab on purpose. Either it is for my research project purposes  or due to my research project purposes. That is why I came up with an idea to do something in Matlab that would not have any application. I wanted to do something what would look nice and be by itself. Later of course I found some interesting application areas of my work but let us start from the beginning.

Łukasz Kocewiak
Photograph by Magdalena Sozańska

I decided to use Image Processing Toolbox because, in my opinion, no matter of what working on images should always give interesting results. Why not to use some of my pictures from the past taken with my old Minolta Dynax 700si. I still should have it somewhere in the basement.

%% Gausian convolution matrix
G= [1  4  7  4 1;
    4 16 26 16 4;
    7 26 41 26 7;
    4 16 26 16 4;
    1  4  7  4 1];
G= G/sum(sum(G));

%% Picture analysis
ExampleUint= imread('Example.jpg');
ExampleDouble= double(ExampleUint);
[m n]= size(ExampleDouble(:,:,1));

%%  Linear Gaussian filter
ExampleFilterUint= imfilter(ExampleUint,G,'conv');
for k=1:10,
    ExampleFilterUint= imfilter(ExampleFilterUint,G,'conv');
    k= k+1; %#ok<FXSET>
ExampleFilterDouble= double(ExampleFilterUint);

Due to the fact that the picture was taken still in analog technology and scanned in any photo lat the quality is not so good. Even in the small picture (a) noise can be clearly seen. I decided to get rig of it by application of a digital filter. A lot of noise and sharp edges that can affect the final result. In order to make the picture more smooth we will use a low-pass filter (b). I decided to use Gaussian filter with the convolution matrix as presented above.

Łukasz Kocewiak 3DThis apprioach can be succesfully used if we would like to evaluate how objects are shaped by light. A good example is in case of body shaping by either attending gym or fitness. If we simply cannot assess if our muscles are sufficiently shaped, just take a picture and do some processing in Matlab. Probably other applications can be found. Please let me know if you have any suggestions by posting a comment.

And an exemplary code from Matlab showing how to display two-dimensional matrix of double precision numbers.

scrsz= get(0,'ScreenSize');
%% Show results
fi1= figure('Name','3D Plot',...
     'Position',[0.1*scrsz(3) 0.1*scrsz(4) 0.35*scrsz(3) 0.5*scrsz(4)]);
ax1= axes('Parent',fi1,'FontName','Verdana','FontSize',10);
view(ax1,[-150 70]),xlim([0 n]),ylim([0 m]),