Technique for the comprehensive characterization of supraharmonic disturbances (9–150 kHz) in the joint time-frequency domain

Sustainable Energy, Grids and Networks 36 (2023) 101181 Available online 27 September 2023 2352-4677/© 2023 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/bync-nd/4.0/). Technique for the comprehensive characterization of supraharmonic disturbances (9–150 kHz) in the joint time-frequency domain Alexander Gallarreta a , * , Igor Fern´ andez a , Jon Gonz´ alez-Ramos a , David de la Vega a , Itziar Angulo b , Amaia Arrinda a a University of the Basque Country (UPV/EHU), Dpt. of Communications Engineering, Bilbao ES-48013, Spain b University of the Basque Country (UPV/EHU), Dpt. of Applied Mathematics, Bilbao ES-48013, Spain ARTICLE INFO Keywords: Electromagnetic interference Measurement techniques Power quality Smart grids Power line communicarions Voltage measurement Supraharmonic impulsive disturbances ABSTRACT The inverter-based devices connected to the low-voltage (LV) grid, such as photovoltaic panels (PVs) and electric vehicle charging stations (EVCSs), are known to generate supraharmonic distortions in CISPR Band A (9–150 kHz) that may interfere power line communications (PLC). Current measurement methods were designed to provide results of the disturbances only in the frequency domain. Thus, impulsive signals and fast amplitude variations of the disturbances are lost in the calculation process. This work addresses this issue by developing a statistical study of the distribution of supraharmonic disturbances in the joint time-frequency domain, which allows the clear differentiation of impulsive signals. The results show that the amplitudes of the normalized measurements can be modeled by the t Location-Scale probability density function (PDF), which can be used to set a threshold for the differentiation of the impulsive disturbances. The study proposes a measurement method and a set of metrics for the complete characterization of the disturbances present in the LV grid, by characterizing the noise and distortions as a whole and also impulsive emissions separately. The main contribution of this work is that the measurement technique characterizes the disturbances in the joint time-frequency domain, and it is coherent with methods below 9 kHz. The proposed method is applied to some typical grid recordings containing different types of noise and emissions. The results demonstrate that the method allows the complete characterization of the supraharmonic distortions present in the LV grid and a proper differentiation of the impulsive emissions in the joint time-frequency domain. 1. Introduction It is widely known that the equipment connected to the low-voltage (LV) grid injects disturbances in CISPR Band A (9–150 kHz). Moreover, inverter-based devices producing or consuming vast amounts of energy, such as photovoltaic panels (PVs) or electric vehicle charging stations (EVCSs), generate high amplitude distortions [1–3]. This equipment has a considerable presence in the grid and is expected to be deployed on a massive scale in the coming years to reduce humans’ carbon footprint, according to the Climate Agreements promoted by the United Nations (UN) or European Union’s policies. In the last decade, several studies have been carried out to classify and parametrize the distortions in this frequency range (commonly named as ‘supraharmonics range’), based on the characteristics of the waveforms [4–7]. Additionally, the impact of these disturbances on the equipment connected to the Low Voltage (LV) grid has been remarked in the literature. The analyzed adverse effects include malfunctions, thermal stress or the production of audible noise in nearby equipment connected to the grid [8–11]. Furthermore, power line communication (PLC) technologies can be interfered by supraharmonic disturbances [12–14], as they use the LV grid as the propagation channel for communications in CISPR Band A, i.e., to send telemetry data generated by the electric smart meters. Moreover, the impulsive disturbances in the LV grid may cause the channel estimation performed by PLC receivers to be incorrect, preventing from equalizing and decoding the frames overlapped with impulsive emissions [15–18]. Although there are numerous PLC technologies available in the market, some of them based on burst transmissions and others on continuous communications, all of them can be interfered with impulsive disturbances. Impulsive distortions in CISPR Band A are common waveforms in the LV grid and are linked to various factors, such as the switched-mode power supplies (SMPSs), the switching frequency of high-power semiconductor components, the pulse width modulation * Correspondence to: Dpt. of Communications Engineering, University of the Basque Country (UPV/EHU), ES-48013 Bilbao, Spain. E-mail address: [email protected] (A. Gallarreta). Contents lists available at ScienceDirect Sustainable Energy, Grids and Networks journal homepage: www.elsevier.com/locate/segan https://doi.org/10.1016/j.segan.2023.101181 Received 30 March 2023; Received in revised form 18 September 2023; Accepted 25 September 2023 Sustainable Energy, Grids and Networks 36 (2023) 101181 2 (PWM), and the short duration and high amplitude transients [18–20]. These distortions are injected by the inverters of the PV panels and the EVCS, which implement SMPS and PWM, at the point of connection (POC) in the LV grid [21–23]. It is at these locations where smart meters are installed to measure the energy consumption and generation of these devices, so that they are located aside the sources of disturbances. Therefore, it is crucial to detect supraharmonic impulsive disturbances, and their characteristics need to be analyzed separately from the other waveforms present in the grid. Impulsive distortions in the LV power grid have been studied for decades and there are numerous contributions in this regard. The majority of the studies are based on the characterization of the temporal behavior of disturbances, by means of metrics in the time domain, such as the period between impulses, the pulse width duration or the amplitude value of impulses [24–26]. Nevertheless, in these studies the disturbances are not differentiated in frequency; therefore, the impulsive supraharmonic distortions are analyzed together with the rest of the impulsive disturbances present in the grid (transients, emissions generated by partial discharges, etc.). The most recent contributions parametrize impulsive disturbances in the time and frequency domains, but separately [27–30], i.e., they do not characterize impulsive distortions in the joint time-frequency domain (on the basis of spectrograms), but considering the time analysis and the frequency results independently. Currently, several measurement methods are available for the characterization of the supraharmonic disturbances (9–150 kHz) propagated through the LV grid. Although there is still no standardized measurement method in this range, the Annex C of IEC 61000–4–30 Ed.3 [31] outlines three temporary proposals: the receiver described in CISPR 16 standard series [32], an extension of the IEC 61000–4–7 Annex B method [33] to frequencies above 9 kHz, and an alternative method specified in the same Annex C of IEC 61000–4–30 [31]. Besides, several research methods have been described in the literature: the Wavelet Approach [34], the Subsampling Approach [35], the RM-A [36], the Light-QP [37] and the Statistical-QP methods [38]. Nonetheless, all these standardized and research methods deal with the characterization of the supraharmonic distortions only in the frequency domain and do not make a differentiated analysis between stable emissions over time or impulsive emissions. Additionally, current methods are based on calculating root-mean-square (RMS) and quasi-peak (QP) values during time intervals equal to or longer than 200 ms, and look for representative values for different aggregation intervals (200 ms, 3 s, 10 min or 1 day) [31]. This implies that impulsive disturbances, with durations of a few milliseconds, are lost in the calculations of averaging and aggregation processes within these intervals. Impulsive disturbances require specific parameters and methods of analysis, different from those seeking representative values for the whole aggregation interval. Thus, in order to simultaneously measure average values while identifying and parameterizing impulsive noises, the metrics and measurement methods for both aspects must be compatible and complementary. The current measurement framework for assessing distortions in the entire supraharmonic region (2–150 kHz) presents a principal drawback since the whole band has no continuity on metrics and amplitude limits. The compatibility levels (CL) regulated in IEC 61000–2–2 [39], which are the maximum amplitude levels for supraharmonic disturbances propagated through the LV in the 2–150 kHz range, are defined with two different and non-comparable metrics depending on the frequency band. In the 2–9 kHz range, the CL are defined in the form of RMS, which must be calculated with the IEC 61000–4–7 Annex B method. The QP values provided by the CISPR 16 method are needed to assess these limits in the 9–150 kHz range. These discontinuities present inconsistencies in the results for the supraharmonic range, as the results below and above 9 kHz cannot be compared. This work proposes a solution for the complete characterization of supraharmonic distortions by first studying the distribution of distortions in the joint time-frequency domain for CISPR Band A. Based on the results obtained in that study, this paper proposes a measurement method for the parametrization of all distortions present in the grid and also for the detection and separate characterization of impulsive disturbances in the joint time-frequency domain, using simple, compatible and complementary metrics to the measurement framework defined for CISPR Band A. Moreover, the proposed method ensures the coherence in result in the entire 2–150 kHz band, since it is based on an adaptation of the IEC 61000–2–4–7 Annex B method for the CISPR Band A. This paper is organized as follows. Firstly, a state of the art of measurement methods for CISPR Band A is presented. Secondly, the objectives of the work are described. Thirdly, the performance of supraharmonic measurement methods for synthetic waveforms with different repetition rates and amplitudes is studied. Afterwards, the methodology of the work is proposed. Then, the results of the distributions obtained for supraharmonic distortions are analyzed, and the conditions for which a disturbance can be considered impulsive are presented. Subsequently, a technique for the complete characterization of supraharmonic disturbances is proposed, in which a separate analysis of the impulsive disturbances is proposed based on the results obtained in the work. After that, the application of the proposed technique to real measurements recorded in the LV grid is shown. Finally, the conclusions of this paper are described. 2. Measurement methods for conducted emissions in the LV grid This section describes the current state of the art of measurement methods to characterize the supraharmonic disturbances in the LV grid. This section describes, firstly, the standard measurement methods for CISPR Band A and, secondly, the research methods available in the literature. 2.1. Standardized measurement methods The measurement methods to characterize the conducted disturbances that are propagated through the LV grid are not completely standardized. As for frequencies up to 9 kHz the standard IEC 61000–4–7 [33] defines metrics and assessment methodologies, for frequencies higher than 9 kHz there is still no standard method. Due to the significant relevance that the analysis of the disturbances in frequencies up to 150 kHz is gaining, the standardization bodies are working in the development of new methods (IEC SC77A WG9) or in the compilation of recent research work (CENELEC TC219 WG11). Therefore, there is a need of defining adequate metrics and developing novel measurement methods in this frequency range. On the watch of the definition of a standard method, Annex C of IEC 61000–4–30 Ed.3 [31] suggests three temporary assessment procedures for measurements in the 9–150 kHz band. Firstly, the method described in CISPR 16 [32] standard series is proposed, whose implementation is based on a super-heterodyne receiver using analog circuitry. This method provides QP values in the frequency domain with a resolution bandwidth of 200 Hz at −6 dB and a high time resolution due to a high degree of overlap in the time windowing of the input values. Nonetheless, this method requires a high computational burden and a large amount of memory resources. Thus, IEC 61000–4–30 Ed.3 standard states that this method is not adapted for on-field surveys; on the contrary, it is a measurement method defined for electromagnetic compatibility measurements, with a laboratory setup defined in the set of CISPR 16 standards that makes its application in field measurements impossible, mainly due to the mandatory use of a line impedance stabilization network (LISN). Secondly, Annex C of IEC 61000–4–30 [31] proposes an alternative method with lighter computational complexity, whose description is provided in the same annex. The method is defined to provide spectral results with a resolution bandwidth of 2 kHz by processing only 8 % of the input signal, i.e., 92 % of the spectral content is omitted in the analysis. These aspects generate a significant loss of A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 3 events in the time domain and a rough frequency resolution, and therefore, this method is not commonly used to assess the amplitude of the disturbances in the LV grid. Thirdly, the use of the procedure described in IEC 61000–4–7 Annex B [33] is proposed to be extended to the 9–150 kHz range, although how to adapt this assessment method to a much wider frequency range is not described in the standard. The IEC 61000–4–7 is a gapless method based on non-overlapped rectangular windows of 200 ms to compute the spectral content of the input data. As in CISPR 16–1–1, the results have a resolution bandwidth of 200 Hz. Although the proposal of measurement methods made in IEC 61000–4–30 is not ‘normative’, none of them addresses the analysis over time. This implies that relevant information, such as impulsive noise and time-varying disturbances, is lost in the assessment process, and therefore, is not shown in the final results. New methods in this band providing not only good frequency resolution, but also the amplitude variation of disturbances in the joint time-frequency domain are needed for a detailed analysis of the supraharmonic disturbances. These characteristics are addressed in the proposed method. 2.2. Research measurement methods Several research methods have been proposed in the scientific literature, such as the Wavelet Approach, the Subsampling Approach, the Light-QP, the RM-A and Statistical-QP methods. The Wavelet Approach [34] analyses the spectral content of the sampled data using the Wavelet Packet Decomposition (WPD), instead of discrete Fourier transforms (DFTs). This procedure filters and downsamples the signal recursively to obtain 740 frequency bins of 200 Hz in the 2–150 kHz frequency range. The amplitude value of the signal for each frequency bin is obtained by computing the root-mean-square of filtered samples in measurement intervals corresponding to ten cycles of the mains (200 ms). The Subsampling Approach [35] is based on the Nyquist-Shannon sampling theorem, which defines that the sampling frequency must be at least twice the maximum frequency of the waveform to be analyzed. This method was designed to characterize supraharmonic disturbances using existing platforms with limited sampling rates. This method implements an analog filter bank to decompose the signal into 10 bands of 15 kHz to analyze the entire 9–150 kHz band. Hence, the maximum sampling frequency required to characterize the conducted distortions with this method is 30 kHz. The output filtered signals are processed with DFTs using rectangular windows of 0.5 ms length. The RM-A method [36] is an adaptation of the IEC 61000–4–7 Annex B for the CISPR Band A (9–150 kHz). This gapless method applies non-overlapped rectangular windows of 20 ms length to compute the DFTs. The resultant spectral values have a resolution bandwidth and a frequency-step-size of 50 Hz, which are grouped by performing an ad-hoc ‘symmetrical grouping’ proposed in this method, to obtain a 200 Hz resolution bandwidth and a 100 Hz frequency-step-size. Then, in this method, the grouped values are aggregated in time to obtain both the maximum and root-mean-square for each 200 ms and 3 s, which are the measurement intervals defined in the IEC 61000–4–30. As the disturbances in the 9–150 kHz range are not generated by distortions of the fundamental, but by the connected devices, the definition of these intervals is not related to the period of the fundamental, but to keep the coherence to the methods below 9 kHz. The Light-QP method [37] is an alternative method to CISPR 16–1–1 that provides QP values with a considerably lower computational burden and memory resources. It is a two-stage method, where the first stage is the RM-A method, which provides root-mean-square (RMS) outputs every 20 ms. These RMS values are processed with a specific digital implementation of a quasi-peak detector, to obtain results similar to the CISPR 16–1–1. Lastly, the Statistical-QP method [38] provides QP spectra, comparable to those obtained with CISPR 16, based on a statistical relationship between the RMS and QP values for CISPR Band A. This method uses the instantaneous RMS outputs of the RM-A method, a statistical analysis based on percentiles to quantify the fluctuations of the disturbances and a simple first-degree polynomial equation to estimate the QP values. Thus, QP values can be calculated offline if the instantaneous RMS outputs of the RM-A method are stored during a measurement. As for the standardized measurement methods, none of the research methods characterize the disturbances in the joint time-frequency domain. None of them either identify and analyze the impulsive emissions in the joint time-frequency domain. As a result, the fast amplitude variations are hidden in the results of these methods, due to the weighting process of the detectors and the temporal aggregation strategies. 3. Objectives This work aims at proposing a set of metrics and assessment methodology for the complete characterization of the supraharmonic disturbances propagated through the LV grid in the joint time-frequency domains within the CISPR Band A. This methodology should characterize the grid distortions in the joint time-frequency domain, as well as to parametrize separately the impulsive distortions in this band that might affect PLC. Three main tasks have been defined in order to achieve this objective. Firstly, the study parametrizes statistically the amplitude variation of the conducted disturbances over time, for each 200 Hz-frequency band in the entire frequency range, and for the measurement intervals defined in power quality (PQ) surveys [31]. This statistical characterization is the basis to define a criterion for distinguishing impulsive disturbances in the time domain. Secondly, a calculation procedure to identify impulsive disturbances is defined, as this type of disturbances may be hidden in the spectral estimation over several seconds. Lastly, the necessary metrics to parametrize all distortions jointly and impulsive distortions separately are proposed. 4. Analysis of current methods All the methods currently available were designed to provide results only in the frequency domain. None of them provide information about the amplitude variations of the supraharmonic distortions over time, but only provide a weighted spectral result in a measurement interval. As a result, different signals with different characteristics over time (repetition frequency, amplitude, etc.) may give similar outputs in the frequency domain. Therefore, it is necessary to define metrics that provide representative results of the variability of the emissions within the measurement interval. Moreover, these metrics should be comparable, or easy to relate, for the most commonly used measurement methods, so that the information provided by these metrics is complementary to the results of the current methods. In order to illustrate this reasoning, the most relevant measurement methods have been applied to three LV grid waveforms. All the waveforms contain a disturbance at 132.2 kHz, in the form of a set of noncontinuous bursts, but with a different number and duration of the bursts: the waveform ‘Grid-1 ′ contains 14 bursts of 76 dB μ V and 54 ms (see Fig. 1), ‘Grid-2 ′ contains 6 bursts of the same duration (208 ms) and amplitude (76 dB μ V) at 132.2 kHz (see Fig. 2); and ‘Grid-3 ′ contains 6 bursts of 80 dB μ V and different duration (between 108 ms and 346 ms) at 132.2 kHz (see Fig. 3). A selection of standardized (CISPR 16–1–1, IEC 61000–4–7) and recently published (Light-QP, RM-A) measurement methods have been applied to the three grid waveforms. Although the variation pattern over the time of the disturbance at 132.2 kHz is noticeably different for the three grid waveforms, the outputs are similar in amplitude for the same detector (RMS or QP), regardless of the method considered in the calculation (see Fig. 1 to Fig. 3): •RMS: 73 dBuV, assessed by IEC 61000–4-7 and RM-A methods. A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 4 •QP: 76 dBuV, assessed by CISPR16-1–1 and Light-QP. This is a representative example that demonstrates that disturbances with different variation patterns of amplitude, duration and repetition rates are evaluated similarly by current standard measurement methods. Novel methods recently published, as RM-A [36] or Light-QP [37], do not evaluate this aspect either. Moreover, the information about the temporal evolution of the disturbances, which can be relevant to evaluate the impact on electronic devices [8] or PLC transmissions [15], or to identify the source of the emission [40], is lost in the assessment procedure of the measurement method during measurement intervals defined in the standard [31]. This fact highlights the need to develop novel techniques and methods to characterize the disturbances in the joint time-frequency domain and to define adequate metrics to parameterize the variability of the distortions and, if it is the case, to identify and characterize the impulsive nature. Fig. 1. Spectrogram and spectrum of waveform ’Grid-1 ′ . Fig. 2. Spectrogram and spectrum of waveform ’Grid-2 ′ . Fig. 3. Spectrogram and spectrum of waveform ’Grid-3 ′ . A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 5 5. Methodology This section describes the methodology followed to parametrize the supraharmonic disturbances in the joint time-frequency domain. Since the measurement method plays a pivotal role in the characterization of conducted emissions, firstly, a methodology to characterize the supraharmonic disturbances in CISPR Band A is presented in this study. Secondly, the procedure for the statistical parametrization of the disturbances in the 9–150 kHz band is described. Lastly, the proposed assessment methodology is evaluated with a dataset of recordings. 5.1. Selection of the measurement method for the statistical characterization This study proposes a processing technique that is consistent with existing methods in the literature, which allows the direct relation between methods and/or the comparison of results. As there is no normative method for the 9–150 kHz range, the existing research and standardized methods have been reviewed to select the most appropriate to be used as a basis for developing a more complete (time-frequency domain) processing technique. Among the existing methods, only those that propose results in the form of RMS are selected, discarding the outputs of the QP, average or peak detectors. The essential reason is that the RMS output is the most widely used metric in PQ, because it can be easily related to intentional emissions in the same frequency band (PLC) [17], to PQ degradation [3] and to physical effects on the devices (thermal stress and malfunction, mainly [11]). In [41], an accuracy analysis of most of the methods mentioned in Section 2 is performed. This concludes that the IEC 61000–4–7 Annex B method allows the most accurate characterization of the distortions in the CISPR Band A with RMS values. Recently, the method RM-A has been presented with the purpose of adapting the IEC 61000–4–7 Annex B method to the CISPR Band A, providing good accuracy and good resolution in time and frequency [36]. Compared to IEC 61000–4–7, in the frequency domain the number of results is twice (from 200 Hz to 100 Hz between adjacent results) and ten times more in the time domain (results every 20 ms, instead of 200 ms), for the same resolution bandwidth (200 Hz). This is considered an optimal configuration for the detection of supraharmonic disturbances. Furthermore, the RM-A method implements the ‘symmetrical grouping’, a novel technique to group the frequency components, which allows the better allocation of high-frequency distortions in the frequency domain. Last, the computational burden of this method is lower with respect to the rest of the methods, due to the use of non-overlapped rectangular windows of 20 ms length. For all these reasons, the RM-A method has been selected as the basis for the development of a novel method that provides results in the joint time-frequency domain, since it is the method that provides better accuracy and coherence with the RMS values from the IEC 61000–4–7 method for frequencies below 9 kHz. According to [36], the RM-A method is divided in three main processing blocks: the frequency analysis, the frequency grouping and the time aggregation (see Fig. 4). The voltage measurements performed in the LV grid are processed, firstly, with the short-time Fourier transform (STFT), using non-overlapped rectangular windows of 20 ms length. Then, the output values (Y ′ C,b ), which have a resolution bandwidth and frequency-step-size of 50 Hz, are grouped by applying the ‘symmetrical grouping’ (see Eq. (1)) to obtain spectral results with a resolution bandwidth of 200 Hz (Y ′ B,b ). Y ′ B,b= 1 2⋅Y2 C,b−100Hz +∑ b+50Hz f=b−50Hz Y2 C,f+1 2⋅Y2 C,b+100Hz √ √ √ √(1) Additionally, this process is repeated every 100 Hz to achieve a frequency-step-size of 100 Hz. Finally, the time series of the grouped spectral values are aggregated by means of the assessment of the RMS and maximum values in 200 ms and 3 s measurement intervals, according to the statement in IEC 61000–4–30 standard for PQ surveys [31]. 5.2. Proposed methodology for the analysis of supraharmonic disturbances in the joint time-frequency domain A measurement interval of 3 s is proposed in this methodology, which corresponds to the second aggregation interval defined by IEC 61000–4–30 standard. This measurement interval has been chosen since the first interval defined in the standard is 200 ms, which is too short to differentiate the impulsive distortions, and the third aggregation period is 10 min, which would imply a high computational cost due to the high number of DFTs that would have to be calculated. The procedure followed to statistically parameterize the distribution of supraharmonic disturbances in CISPR Band A is as follows: •Data normalization: So as to compare all the available frequency bands obtained with the RM-A method, the values should be normalized. This normalization is performed by computing the median and subtracting this value (in dB μ V) to all the data per each frequency band. Thus, the resulting normalized values are expressed in dB. •Obtaining the optimal parameters for each distribution: In order to obtain the distribution of the normalized data, the parameters of the probability density functions (PDF) that provide the maximum likelihood estimation (MLE), which is an estimation method to compute the parameters of a PDF that better fit the input sample set, are calculated. •Selection of the PDF that best fits the input data: The goodness of fit of different PDFs are compared by computing the Chi-Squared value, which is used to evaluate the hypotheses about the distribution of the observation data. The selected PDF in this study would be the one that provides the highest MLE for its distribution and the lowest value in the Chi-Squared value. Fig. 4. Schematic overview of the RM-A measurement method. A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 6 6. Statistical analysis of the disturbances in the joint timefrequency domain This section describes the analysis performed to parametrize statistically conducted disturbances in the joint time-frequency domain. For that purpose, the PDF that better fits the amplitude variation of the RMS values of these disturbances over time and for each 200 Hz frequency band in the entire frequency range is computed. 6.1. Dataset of recordings measured in the LV grid The statistical analysis is computed with voltage recordings measured in several locations of the LV grid. The dataset is composed of 71 signals measured in the POC of smart meters (49 recordings), EVCS (16 recordings) and photovoltaic inverters (6 recordings). All recordings used in this work were taken in 3 measurement campaigns in Spain. All POC were located between several tens and a few hundred meters from the transformers. Thus, most of the measured disturbances were generated by the equipment or the house connected to the POC of the LV grid. The measurements taken at the POCs of the smart meters were measured in the cabins where these devices are located inside apartment buildings. These measurements were taken in three geographic areas linked to locations dependent on LV transformers; different grid topologies and homes density can be found in these three geographic areas. These measurements were carried out in urban and rural areas, of which 33 recordings correspond to the first scenario, and the remaining 16 to the second one. The EVCSs used in this work were three-phase AC devices with a maximum power of 22 kW. These measurements were taken in 4 locations using 4 different EVs, whose battery level were between 75 % and 90 %. In addition, all PV recordings were measured in a unique installation with 3 connection configurations of 4 panels, using an inverter of 18 kW and different loads. All the measurements were recorded with the same measurement system, which is described in [42]. A voltage probe is used to acquire the signal in the field, which applies a galvanic isolation and a band-pass filter to protect the measurement equipment against overvoltage transients and to filter the emissions outside the 10–500 kHz band [43]. The output of the probe is sampled with a sampling rate of 8.92 MHz and a 16 bit/sample resolution by a digital oscilloscope and recorded in a laptop. 6.2. Statistical analysis The basis of the statistical analysis is the RMS values of the spectrograms, which represent a detailed assessment of the noise and disturbances in the joint time-frequency domain for the measurement interval and for the whole frequency range, as it can be observed in Fig. 1, Fig. 2 and Fig. 3. As it is justified in Section 5.1, the method selected to assess the spectrogram is the RM-A method (Y ′ B,b in Eq. (1) and Fig. 4), which provides RMS outputs in the 9–150 kHz band, with a spectral value every 100 Hz, in 3 s measurement intervals, with a time granularity of 20 ms. This corresponds to a matrix of 1409 frequency bands and 150 temporal values; therefore, 211,350 values are generated for each signal, which implies a total of 15,005,850 values for the entire dataset of 71 recordings. Fig. 5 shows the distribution of the 15,005,850 normalized values of the entire dataset. The normalized values have been modeled with 23 distribution functions: Beta, Binomial, Birnbaum-Saunders, Burr, Exponential, Extreme Value, Gamma, Generalized Extreme Value, Generalized Pareto, Half-normal, Inverse Gaussian, Logistic, Loglogistic, Lognormal, Nakagami, Negative Binomial, Normal, Poisson, Rayleigh, Rician, Stable, t Location-Scale and Weibull. Nevertheless, only 5 out of the 23 distributions have provided valid results, since some of the distribution functions do not converge and others are not able to model negative values. In the distributions providing valid results, the configurations that best fit the normalized data have been obtained, which has been numerically verified, as these configurations provide the highest value of MLE (see Section 5.2). The resultant distribution functions are shown in Fig. 6. The goodness of fit of different PDFs are compared by means of the Chi-Squared value, as it reports about the suitability of the hypotheses for the observation data. The lowest Chi-Squared value represents the most appropriate distribution function for the dataset. Table 1 contains the values of the parameters of each PDF that give the highest MLE and the Chi-Squared results for that configuration. The results in Fig. 6 and Table 1 show that the t Location-Scale is the PDF that best represents the distribution of the normalized input data, since it provides the lowest Chi-Squared value. The t Location-Scale distribution is defined with Eq. (2) [44], whose cumulative distribution function (CDF) is represented in Fig. 7. Fig. 5. Distribution of the normalized spectral values. Fig. 6. Distribution of the normalized dataset and the PDFs with the configuration with highest MLE. A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 7 f(x; μ , σ , ν ) = Γ( υ +1 2) σ  υπ √Γ( υ +1 2)(1+1 υ (x− μ σ )2)− ( υ +1 2)(2) Where x is the input, μ is the location variable, σ is the scale parameter and ν is the shape variable. As a result, the CDF shown in Fig. 7 is the cumulative distribution that best represents the normalized disturbances in the LV power grid. The impulsive disturbances are linked to the highest values of this CDF (close to the 100 %), as they are very infrequent occurrences of amplitude considerably greater than the median value, corresponding to 50 % of the CDF. For that reason, the impulsive disturbances are represented by the upper tail of the distribution, for a small percentage of amplitudes several dBs greater than the median. This CDF can be used as a basis to identify the impulsive disturbances in the grid measurements. For this purpose, three threshold values have been defined to differentiate the occasional high-amplitude impulsive disturbance from the vast majority of the spectral components of stationary waveforms: 90th, 95th and 99th percentiles of the CDF of the normalized values. Table 2 shows the amplitude difference of the percentiles 90th, 95th and 99th with respect to the median value. The results obtained in Table 2 show a considerable difference for close percentiles, since these values are at the upper tail of the distribution; therefore, the selection of the threshold for identifying impulsive distortions depends on how stringent is the criterion applied to detect the supraharmonic disturbances more precisely. The results of applying each of these threshold values is analyzed in the following paragraphs. 7. Proposed metrics and procedure to identify and parametrize impulsive grid disturbances The results obtained in the previous sections show that: •The background noise and conducted disturbances present in the power grid can be statistically modeled as a t Location-Scale PDF. •The CDF of this PDF can be used as a basis to differentiate impulsive disturbances in the 9–150 kHz band, as they correspond to the upper tale of the CDF (a small percentage of sporadic occurrences with an amplitude several dBs higher than the median value of the PDF). •A threshold in the upper tail of the CDF is proposed as a statistical procedure to differentiate the impulse responses in conducted disturbances in the grid. In this work, three percentiles of the CDF have been proposed as a first approach (90th, 95th and 99th). •Once the impulsive disturbances are identified, a set of specific metrics can be used to characterize their relevance, duration and repetition rate. •The RM-A method has been identified as the most appropriate basis for developing a method in the joint time-frequency domain, with metrics defined to characterize the supraharmonic disturbances and identify and evaluate the impulsive disturbances separately. Moreover, the proposed method maintains the coherence with the RMS values obtained for frequencies below 9 kHz Based on these results, a procedure for the complete characterization of emission and for detection and parametrization the impulsive distortions in the CISPR Band A is proposed and described in this section. 7.1. Proposed method to identify and parametrize impulsive emissions in the power grid The procedure for identifying and characterizing the supraharmonic disturbances consists of 3 main blocks (see Fig. 8). Firstly, the RM-A method, described in Section 5.1 is applied to perform the spectral analysis of the raw data. Secondly, the detection stage distinguishes the impulsive distortions from the non-fluctuating emissions over the 3 s measurement interval. Finally, the characterization stage characterizes the amplitude and temporal behavior of the impulsive emissions in each frequency band. Regarding the first stage of the procedure, the spectrogram of the measured voltage (Y ′ B,b ) is obtained by the RM-A method, with a resolution bandwidth of 200 Hz, a 100 Hz frequency-step-size and a time step of 20 ms (see Fig. 9). In the detection stage, the disturbances fluctuating over time are calculated. For this purpose, the median values of the spectra are calculated for each frequency band and for 3 s measurement intervals (U_median 3 s ) and subtracted to the recording. The metric U_median 3 s represents the median value of the disturbances that are present at each frequency band, within the measurement interval of 3 s. The results of the subtraction are the fast variations of the conducted distortions and emissions around the median value. The impulsive disturbances (U_impulsive 3 s ) are calculated by applying the threshold values defined in Table 2, related to a specific percentile of the CDF of the t Location-Scale distribution (3.27 dB, 4.90 dB and 10.55 dB, which correspond to 90th, 95th and 99th). Fig. 10 shows the results of the U_impulsive 3 s for the three thresholds; the spectral components below the threshold are displayed in dark blue (see Fig. 10): Table 1 Parameters of the highest likelihood PDF and the corresponding Chi-Squared value. PDF Parameters Chi-squared value Extreme Value μ =3.3570 σ =10.0686 – 2.1255 ×10 82 Generalized Extreme Value k= − 0.0204 σ =3.3218 μ = − 1.2343 2.4924 ×10 31 Logistic μ =0.0032 σ =1.7610 – 1.2081 ×10 7 Normal μ =0.3806 σ =4.5457 – 2.7062 ×10 25 t LocationScale μ = − 0.0371 σ =1.8777 ν =2.3611 8.0053 Fig. 7. CDF of the t Location-Scale distribution for the normalized input data. Table 2 Normalized values of the t Location-Scale distribution for 3 cumulative probabilities. Cumulative probability (%) Normalized values with respect the median (dB) 90 3.27 95 4.90 99 10.55 A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 8 •The results for the threshold of 3.27 dB (90th of the CDF) show that it does not correctly select the impulsive distortions, since numerous frequency components of the background noise of the LV grid remain as impulsive disturbances. An example of this situation is highlighted with a white oval in the top graph of Fig. 10. •The threshold of 4.90 dB (95th of the CDF) shows a more discriminating selection of the impulsive distortions, but a few frequency components corresponding to background noise still remain. Representative examples are identified with white ovals and a white arrow in the central plot of Fig. 10. •The results for the threshold of 10.55 dB (99th of the CDF) are composed only of impulsive waveforms and all the background noise components are properly filtered out (see the bottom graph in Fig. 10). This plot does not show any spurious frequency components that may correspond to the intrinsic fluctuations of the background electromagnetic noise, which are not impulsive disturbances at all. Thus, this graph only contains emissions that have clearly an impulsive behavior. Therefore, the threshold of 10.55 dB (99th of the CDF) allows the discerning differentiation of impulsive disturbances, which has been tested with the 71 recordings of the data set. Thus, it is proposed to use only the threshold of 10.55 dB in order to correctly detect impulsive emissions. The proposed configuration has been selected to avoid the occurrence of false positives, as this procedure ensures that the outputs of the ‘detection stage’ correspond to impulsive disturbances. This threshold avoids including in the ‘characterization stage’ any remainder of background noise, as demonstrated in Fig. 10. A potential drawback of this criterion is that it may generate false negatives in some Fig. 8. Schematic overview of the procedure to detect and characterize the impulsive disturbances. Fig. 9. Spectrogram of a recording measured in the LV grid obtained with RMA method (Y ′ B,b ). Fig. 10. Waveforms classified as impulsive disturbances (U_impulsive 3 s ) for the thresholds defined in Table 2. A. Gallarreta et al. Sustainable Energy, Grids and Networks 36 (2023) 101181 9 circumstances since the selected amplitude threshold is very high. Nevertheless, this is not considered a disadvantage of the method, as the purpose is the detection of the highest amplitude disturbances. Using the median value of each frequency bin and the threshold corresponding to the 99th percentile of CDF is a stable basis for the detection, ensuring that only the disturbances that fulfill this condition have impulsive behavior. Finally, in the characterization stage, a set of metrics are calculated in order to characterize the impulsive waveforms detected in the previous stage. These metrics are defined in the next section. 7.2. Metrics for the characterization of the impulsive disturbances A set of metrics have been defined, in order to characterize the impulsive noise, in terms of amplitude, durations and repetition during the observation interval. The approach of this characterization is based on the selection of simple, but representative parameters, which give a clear and direct view of the impulsive emissions. Regarding the amplitude characterization, three metrics have been defined to parametrize the impulsive distortions (see Fig. 11): •Max U_impulsive 3 s : maximum amplitude of the impulsive disturbances in a 3 s measurement interval. •RMS U_impulsive 3 s : RMS amplitude of the impulsive distortions over the 3 s observation period. •Min U_impulsive 3 s : minimum amplitude of the detected impulsive waveforms on the 3 s interval. Furthermore, three additional metrics have been defined to characterize the disturbances propagated through the LV grid (see Fig. 11): •U_median 3 s : median value of each frequency band in the 3 s interval. This metric provides an indication of the mean behavior of distortions, since it represents the 50th percentile of the amplitudes. •RM-A: U_RMS 3 s (RMS 200 ms ): RMS values in the 3 s interval of the RMS values aggregated every 200 ms, obtained with the RM-A method. •Light-QP: QP: QP values calculated with the Light-QP measurement method. For the characterization of the emissions for each frequency band over the measurement interval, three metrics are calculated (see Fig. 12): •Num. events: number of occurrences of impulsive disturbances within the 3 s measurement interval. •Total dur.: total time in which the impulsive disturbances are present at each 200 Hz frequency band. •Mean dur.: mean duration of the disturbances at each frequency band is calculated by dividing the total duration of the impulsive disturbances by the number of occurrences identified in that interval. •Mean dur. between events: average duration of the time spacing between supraharmonic impulsive disturbances. This metric only shows results for impulsive disturbances whose total duration is equal to or greater than 100 ms, in order to analyze the temporal spacing of the emissions with the greatest presence in the grid. 8. Application of the proposed method This section shows the application of the measurement method for characterizing the supraharmonic emissions in the grid in the joint timefrequency domain, with a specific analysis to identify and parametrize the impulsive distortions. Two measurements in the LV grid where PLC transmissions are present and impulsive distortions can cause interference in PLC transmissions are used for this purpose. The first recording (’Rec-1 ′ ) contains PLC transmissions (highlighted with black arrows in Fig. 13), in the form of bursts of a few milliseconds in the frequency range 42–89 kHz (delimitated with two black lines in Fig. 13). Moreover, numerous impulsive distortions occupying almost the entire CISPR Band A can be seen in the spectrogram of Fig. 13 (indicated with a white oval). These distortions show high amplitude and fast time variability within the frequency band designated for PLC, and therefore, they may interfere the reception of the PLC signals. ‘Rec-1 ′ has been analyzed with the technique proposed in this work, and the results of the joint time-frequency characterization are shown in Fig. 14 and Fig. 15. The impulsive distortions contained in this recording have been differentiated and assessed separately by specific metrics: •The metric ‘RMS U_impulsive 3 s ’, which represents the RMS values only for the impulsive distortions, shows a ripple in amplitude between 40 kHz and 85 kHz; this implies that the impulsive distortions are of higher amplitude than the PLC transmissions, which usually show a flat amplitude spectrum. •The metrics in the joint time-frequency domain clearly identify the impulsive distortions in the total number of events, where a higher number of impulsive events are detected in the same frequencies where the ripple in the RMS values is observed. •Moreover, the different values in the total duration of events or in the mean duration between events reinforce the fact of the presence of the impulsive emissions alternating the PLC transmissions (see Fig. 15). The absence of the analysis in the joint time-frequency domain, this is, the use of only the RMS or QP spectra of the whole recording in Fig. 15, does not reflect the impulsive interfering emissions at all, as they are hidden by the assessment process. The recording ’Rec-2 ′ also contains PLC transmissions in the 42–89 kHz range (delimitated with two black lines Fig. 16 and Fig. 17) and several impulsive distortions. Around 130 kHz, 6 impulsive emissions in the form of bursts of similar duration can be appreciated (see white rectangles in Fig. 16). The amplitude characterization and the time-frequency behavior of the recordings are shown in Fig. 17 and Fig. 18, respectively, where black boxes have been included to highlight the impulsive distortions. Results in Fig. 17 do not reflect if the emissions in 130 kHz are constant over the time or not, as its variation over the time is lost. On the Fig. 11. Amplitude characterization of the disturbances and the metrics for impulsive distortions. A. Gallarreta et al.