A Method for Real-Time Prediction of the Probability of Voltage Sag Duration Based on Harmonic Footprint

In this paper, a novel method for real-time prediction of voltage sag duration is proposed. It is based on the recently introduced new characteristic of voltage sag, named harmonic footprint, and is formulated using a logistic regression model. The concept is mathematically formulated and statistically analyzed using an extensive set of real grid measurement data which are recorded in distribution grids. Furthermore, the proposed method is applied as a part of an advanced grid-tie converter control. It is included in previously developed methods for fast sag detection and magnitude of voltage sag prediction. The algorithm is applied to the control of grid-tie converters used in distributed generators and tested with real/grid measurement data in the IEEE 13-bus test grid by simulations and in the IEEE 33-bus test grid using a hardware-in-the-loop (HIL) microgrid laboratory testbed. It is shown that this method can prevent unnecessary tripping of distributed generators (DG) and improve low-voltage ride-through (LVRT) support. In addition, the model has the potential to be applied to a wide range of devices or algorithms for the protection, monitoring, and control systems of distribution grids.


I. INTRODUCTION
Modern distribution networks are characterized by increased penetration of distributed energy resources (DERs), that is, distributed generators (DGs) based mainly on renewable energy sources (RESs). The Distribution System Operator (DSO) or Distribution Network Operator (DNO) has a major role in managing energy from the generation sources to the consumers, in ensuring reliability and efficiency in the operation of systems that have DERs, and in providing the quality of supply (specified electric energy quality or power quality (PQ) parameters and quality standards) for servicing system users [1].
Voltage sags or dips are among the most important PQ parameters in contemporary distribution grids. They are a consequence of faults in the power grid, high-power The associate editor coordinating the review of this manuscript and approving it for publication was Arash Asrari . induction motors starting, utility transformer energizing, etc. Their effects can be observed over a wide area, affecting a large number of loads and connected distributed energy resources [2], [3]. Voltage sags are represented by a momentary decrease (10% to 90%) in root-mean-square (RMS) value of the AC voltage (at one, two, or all lines) at the power frequency of duration of 0.5 cycles to 1 min [4], or defined as a sudden reduction in the voltage at a point in the electrical system, with a duration between half a cycle to a few seconds [5]. Negative effects are well-known and may be followed by significant financial losses [6], [7]. They range from a disturbance in the supply of PCs or similar low-power electronic devices to problems in the operation of computers or data centers, industrial drive tripping, or discontinuation of DER generation. Hence, they can both be seen on the demand and generation sides.
On the demand side, to sustain operation and prevent load tripping-off, voltage-tolerant or ride-through curves VOLUME 10, 2022 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ were set for computer (microprocessor-based) equipment by ITIC/CBEMA [8], [9]. However, fulfilling these requirements is insufficient to prevent negative effects. To mitigate them, different types of low-power uninterruptible power supplies (UPSs) [10], or utilization of PQ conditioning devices, such as Ferroresonant transformers [11], Dynamic Voltage Restorers (DVRs) [12], Unified Power Quality Conditioners (UPQCs) [13], Universal Power Quality Conditioning Systems (UPQCSs) [14], distribution static compensators (DSTATCOM) [15], or others have been proposed. It is important to emphasize that for their proper operation voltage sag detection, voltage reference generation (synchronization) and adequate control are essential [16]. Different methods have been proposed for the aforementioned devices [16], [17]. However, there is room for improvement, especially regarding the detection speed, full characterization, and control robustness. The DNO role is connected to energy management, that is, proper operation of the generation side and energy transfer to the users, so different strategies are used (centralized, decentralized, decentralized multi-agent-based, etc.), especially in the case of faults [18], [19]. If there are instantaneous, momentary, or temporary PQ problems (such as voltage sags), the DNO should react promptly, as it will affect the generation ability of the DER. Besides protection of DERs (voltage sags may provoke current overload, control circuit miss operation due to loss of synchronization, an increase in DC link voltage, DC power oscillation, decrease in the lifetime of power switching devices and DC link capacitors, etc. [20], [21]), the DNO may have some specific or even opposite requirements for operation of DERs during a fault. They are given in national grid code rules and respective standards to enable low-voltage ride-through (LVRT). They are stipulating that a DER should stay connected during a voltage sag for a predefined short time and support the grid by providing a certain amount of reactive power (Q) and active power (P) depending on the voltage sag depth and duration [22], [23]. This task requires momentary adaptation of the DER control system and the setting of a new operation point [24]- [27]. As swift voltage sag detection and characterization are required, improvements to existing methods are needed, again.
Voltage sags are commonly characterized by magnitude and duration which can be easily extracted from standard RMS measurements [28]. Other interesting features include imbalance, phase-angle jump, etc [3]. For the interconnection of modern DERs to the grid, detection of voltage sags is vital, as is the determination of their magnitude (the lowest voltage in all phases), unbalance, and duration (the longest duration measured across all channels) [24], [25].
In previous studies, the authors introduced a harmonic footprint method for voltage sag characterization. It demonstrated superior performance for voltage sag detection, magnitude prediction, and imbalance [40], [47]. It was shown, that the voltage sag detection time could be shortened to less than 1 ms [40]. Further research revealed that additional acceleration can be achieved if the harmonic footprint is used for the magnitude of the voltage sag (MoVS) prediction. It was shown that MoVS can be predicted with 95% accuracy and an average time gain for the early determination of 292 ms [47]. However, the voltage sag duration was not treated.
The voltage sag duration depends on the fault-clearing time provided by electrical protection in the power system (current-limiting fuses, expulsion fuses, distance relays, differential relays, overcurrent relays, etc.) or on the acceleration time of a large AC load or transformer energizing time. Usually, the fault clearing time is in the span of 10 ms to 2 s [48], so sag duration ranges from half a cycle to 1 minute [4], [5]. However, in practical meshed networks with differently tuned protective relays, in cases of application of fuse saving strategy and in grids where there are several connected DERs, the duration of the voltage sag at a DER interconnection bus may not be straightforward. Therefore, introducing a method for predicting voltage sag duration may be useful.
In this paper, the authors present a probabilistic method for the real-time prediction of voltage sag duration based on the harmonic footprint. The main idea is to improve the DER operation by anticipating the sag duration during the transient stage of a fault. The harmonic footprint (selected set of harmonics) was used as the main sag feature to achieve this solution. To the best of our knowledge, no real-time method has yet been published. The establishing of the mathematical model and verification using a set of voltage sags data obtained by real-grid measurements are presented. Such a solution is implemented for DER control circuit improvement and tested in a modified IEEE 13-bus test grid. Further validation was performed using a hardware-in-the-loop (HIL) microgrid system with the IEEE 33-bus test grid, where the method was implemented in a wind turbine grid-tied inverter. It is shown, that the method can prevent unnecessary tripping of DERs and improve LVRT support. However, there are other possible applications such as in the control units of DVRs, UPQCs, UPQCSs, DSTATCOMs, etc., but they will be addressed in future work.
The main contributions of this paper are as follows: • A method for the real-time prediction of probability (RTPP) of voltage sag duration based on harmonic footprint has been proposed; • The proposed method was verified by testing using real grid measurement data; • It was applied in a DER control unit; • The usefulness of the proposed method was tested using simulations and a laboratory Microgrid HIL testbed on IEEE 13-bus and IEEE 33-bus test grids, respectively, and • A comparison of a DER operation with and without proposed method implementation was performed and benefits are noted. The paper is organized as follows: the background and literature overview are presented in Section II. The proposed method for RTPP of voltage sag duration, with mathematical analysis, examination of model effectiveness, and discussion are given in Section III. The RTTP algorithm with model implementation, including early detection and magnitude of voltage sag prediction, is described in Section IV. Applications of the proposed method for the grid-tied inverter operation, and application in the IEEE 33-bus Microgrid HIL setup, and for protection circuit are presented in Section V. Testing of the proposed method using real-grid measurement results by the MATLAB/Simulink model of the IEEE 13-bus test grid and using an advanced Microgrid HIL-based realtime laboratory testbed with an implemented IEEE 33-bus test grid, as well as a discussion of the presented results, are presented in Section VI. The rest of the paper is composed of conclusions, acknowledgments, and references.

II. BACKGROUND AND LITERATURE OVERVIEW A. VOLTAGE SAGS DURATION CLASSIFICATION IN STANDARDS
Classification of voltage sags according to the magnitudeduration plane, as recommended by IEEE 1159 [4], divides sags into three groups: instantaneous with a duration between 8.3 ms and 0.5 s (or 10 ms and 0.6 s in case of 50 Hz AC networks), momentary with a duration between 0.5 s and 3 s (or 0.6 s and 3 s in case of 50 Hz AC networks) and temporary sags with a duration between 3 s and 60 s.
If LVRT requirements are considered, the sag is either instantaneous/momentary or temporary. According to IEEE Std. 1547-2018 [23], and following LVRT recommendations (grid codes) [22], a protective device will trip only if the sag is deep (less than 55% of rated voltage and with a duration of above 160 ms). In other cases, it may ride through or trip, or it has permissive or mandatory operation capability. Therefore, there may be two possible states: stay connected (in cases of instantaneous/momentary sags of duration less than 3 s), or disconnection is allowed (in cases of temporary sags of duration above 3 s). Similar actions should be taken if the German grid code or some other grid code is applied, i.e., if a sag falls above or below the LVRT line. After 3 s, the DG is not obligated to stay connected to the grid, regardless of the MoVS, both in the IEEE recommendations and in the grid code requirements. In addition, because of the high stress on the inverter components, disconnection is likely justified after 3 s, even if LVRT may continue.

B. HARMONIC FOOTPRINT
The novel real-time method proposed in this study is based on the harmonic footprint (HF), a recently introduced voltage sag indicator [40]. The HF is a voltage sag pattern obtained by tracking the voltage harmonic distortion represented by the sum of the 2 nd , the 3 rd , the 5 th , and the 7 th harmonics (HDU2357) during the transition phase of the voltage sag. It usually appears twice during a fault: at the beginning (starting) phase and at the end of the recovery phase of the sag.
The HF has four main characteristics: 1) initialization or starting phase (exponential rise) suitable for detection, 2) maximum value (MHDU), 3) width, and 4) duration. They carry information about the voltage sag's beginning (starting) and magnitude, whereas its duration can be predicted as the probability of developing into a momentary or temporary one, which will be presented further in this paper.

C. LITERATURE OVERVIEW OF VOLTAGE SAGS PREDICTION METHODS
Usually, the determination of the duration of voltage sags is a part of the methods developed for voltage sag prediction. These methods predict the frequency of a voltage sag occurrence during a given time interval and its typical characteristics (magnitude and duration). Computer simulations and long-term measurement data analyses were performed.
Stochastic assessment is the most frequently applied in simulations and includes the critical distance, fault positions, and Monte Carlo (MC) methods. It combines deterministic results with stochastic data to produce a probabilistic assessment of the voltage sags at a given bus in a system [49]- [53]. In [49], a short overview of prediction methods, from essential ones (voltage divider and critical distance) to advanced stochastic methods, is presented. Detailed stochastic models are required to achieve good results and to consider the second-order events. In [50], it was concluded that the method of the critical distances yields acceptable results and is far more accurate than the method of fault positions, but it requires high computational effort. In [51], a state-of-theart offline method based on the stochastic prediction of the duration and magnitude of the voltage sag was presented. The results show a higher accuracy of the presented method in comparison to the traditional techniques. However, this method is considerably more complex. In [52], the application of stochastic-based methods for network planning was addressed. It successfully predicts proper fault parameters, whereas the entire method for grid planning shows ''adequate'' results. MC simulations were used to predict the voltage sags [53]. It is shown that the MC simulation provides a more comprehensive characterization than the fault position (location) method.
Another approach is to use the recorded voltage sag data obtained from long-term measurements in a real grid. After statistical analysis, the sags' duration and magnitude distributions, and the probability distribution functions were obtained [54]- [57]. Recently, such an approach has been used to develop a stochastic model for forecasting voltage sags in a real grid [58]- [61]. In [58], the stochastic characterization of voltage sag data obtained from long-term monitoring in a Portuguese grid is presented. It is shown that the occurrence of voltage sags as a Poisson process is not valid if voltage VOLUME 10, 2022 sags clusters (grouped sags) are included. In the opposite case, a gamma distribution fitting was confirmed. In [59], some ideas for voltage sags forecast using the parameter time to the next event as the probabilistic variable are given. Further elaboration is provided in [60]. In this paper, a gamma distribution is proposed to properly model voltage sags based on the forecast. High accuracy was achieved, i.e., the forecast presented an error slightly greater than 10%. The method was further improved in [61] by introducing special indices of intermittence to detect and remove grouped voltage sags (clusters). In this case, rare voltage sags can be predicted. The results showed that the rare voltage sag forecast had acceptable errors and was not prohibitive.
However, none of the previously described methods provide predictions in real-time, and they require either knowledge of the grid parameters or long-term measurements. Therefore, they are not suitable for application in the control systems of modern DERs and of the devices for PQ improvement.

III. PROBABILISTIC PREDICTION OF VOLTAGE SAGS DURATION A. REAL-GRID MEASUREMENT RESULTS ANALYSIS
To define the method for real-time (online) prediction of the probability of the voltage sag duration (RTPP), the authors analyzed various measurement results presented in the literature [3], [40], [54]- [56]. As the most convenient, the available set of real-grid measurement data, measured at different voltage levels and in different distribution grids, containing 680 recordings, out of which 95 were voltage sags, was selected for further analysis. A detailed description of these voltage sags is given in [40], [47]. The main characteristics of these sags and MHDU are presented in Appendix (Tables 6 and 7), while the magnitude-duration graph is shown in Fig. 1. It can be observed that most of the sags end in the period up to 3 s and that they can be regarded as momentary sags (83 sags / 87%). Very few of them have a longer duration (above 3 s), i.e., temporary sags (12 sags / 13%), which follows the previous conclusion regarding the distribution of results. This observation motivated the authors to perform further analyses using the HF to enable real-time sag duration prediction.

B. HARMONIC FOOTPRINT ANALYSIS
All available measurement records were represented using HF characteristics, but only MHDU was considered. The descriptive statistic of HF regarding MHDU for a sample of size N=95 is presented in Table 1. Momentary and temporary sags are distinguished and the mean value of the MHDU and 95% Confidence Interval (CI) of the MHDU were calculated. It can be observed that momentary and temporary sags can be distinguished and that a probabilistic model may be formed. Therefore, we propose the following hypothesis: Hypothesis: Based on the maximum value of the harmonic footprint (MHDU) during a transient, it is possible to  determine the probability that a voltage sag will develop as temporary or momentary.

C. MATHEMATICAL ANALYSIS -LOGISTIC REGRESSION MODEL
To treat the hypothesis in this research, among well-known traditional techniques (such as ordinary least squares regression or linear discriminant function analysis), the logistic regression model was chosen. Generally, logistic regression is widely used in modern statistics to test hypotheses where the outcome is a categorical variable (in this case, dichotomous), while predictor variables (one or more) can be categorical or continuous [62]- [65]. In practice, many problems require analysis and predictions based on dichotomous variables [66], [67]. Therefore, the probability that a voltage sag is temporary or not can be computed using the logistic regression model.
Let us indicate with X a random variable, which represents MHDU for a given sag with a set of values from 0% to 100%, and with Y a dichotomous variable, which takes a value of 1 for temporary sags (≥3 s) and 0 for momentary sags (<3 s). Then, a probability function of happening Y, for a specific value X = x can be defined as Logistic regression is based on the logit function of one or more dependent variables. The logit function of f (X ) is the natural logarithm (ln) of odds, which represents the ratio of probability f (X ) = P (when Y = 1 for a given value of X ) to probability 1 − f (X ) = P (when Y = 0 for a given value of X ), i.e., The simple logistic model (with one independent variable X) has the form where α is the constant of the model, and β is the regression coefficient. Using equation (3) and taking the antilog of both sides, it is obtained Therefore, the standard logistic function with one independent variable X is defined as where coefficients α and β are estimated from given data.
To calculate the coefficients and further test the performances of the method, the available set of measured data was divided into two groups: main (83 sags) and control (12 sags). Using statistical software Minitab and 95% confidence level (CL), and based on 83 measurements of the main group, the obtained coefficients for the model (presented by (5)) are α = −9.15 and β = 0.23202. The coefficient α is constant and determines the horizontal position of the logistic regression curve. The slope of the curve is specified by coefficient β. As β > 0, the probability that the sag is temporary rises with the increase in the maximum value of the HF (MHDU).
The graph of the logistic curve, i.e., the probability function of the MHDU, is determined by and presented in Fig. 2. The fits and diagnostics for all observations are provided in the Appendix ( Table 8). As it can be seen, the graph is a sigmoidal (or S-shaped) curve, with probability values tending to 0 for X < 30 and 1 for X > 50. Thus, this hypothesis was confirmed. However, fine-tuning is required.

D. FINE-TUNING OF THE RTPP
The calculated value f (X = x) is the probability that the sag with MHDU=x will develop as temporary, whereas 1-f (X ) is the probability that it will be momentary. To obtain fine-tuning of the border (uncertain) zone for the RTPP (Fig. 2), several probabilities were tested in the range of 30% < X < 50%.
As f (X = 30%) = 0.1007, for MHDU values less than 30%, it can be stated that the considered sag will develop as momentary with a probability greater than 90%. As f (X = 50%) = 0.9206, for MHDU values greater than 50%, it can be stated that the considered sag will develop as temporary with a probability greater than 92%. Moreover, as f (X = 52.5%) = 0.9540, i.e., for MHDU values greater than 53%, the sag will develop as temporary with a probability greater than 95%.
As f (X = 37.7%) = 0.4006, f (X = 39.5%) = 0.5073, and f (X = 41.2%) = 0.6009, the values of the MHDU around 40% are border values, so the probabilities of sag duration are approximately the same and the prediction cannot be exactly obtained. Therefore, the values of the MHDU between 37.7% and 41.2% will be treated as an ''uncertain zone,'' and the probability function curve can be updated (see Fig. 2).
Finally, the effectiveness of the model must be examined.

E. EXAMINING THE MODEL EFFECTIVENESS
There are several ways to examine the effectiveness of a model [68]. In this research, the effectiveness of the proposed model is examined with: 1) significance of predictor coefficient, 2) values of R 2 , and 3) measure of goodness-of-fit of the model. The standard Wald test was used to examine the significance of the predictor coefficients [69]. The proposed null hypothesis for this test is that some coefficients are equal to zero (i.e., no effect in the model). In Table 2, besides the values of constants, the following are presented: 95% confidence interval, Z -values (a statistical measurement that represents a value related to the mean of a group of values), and p-values (the probability of finding the observed, or more extreme results, under the null hypothesis). As the p-values for α (constant) and β (coefficient of X ) were less than 0.05, both coefficients were significant in the model. Thus, there is a linear relationship between the logit(f(X)) and X .
An additional measure for evaluating the performance of the proposed model performs is the coefficient of determination R 2 (pseudo-R-squared   be R 2 CS = 0.6881 and R 2 N = 0.6742, respectively. Because these values were greater than 0.5, the model effectively matched the data. Finally, the goodness-of-fit of the model was verified using the Pearson and Hosmer-Lemeshow test [62], resulting in p-values of 1 and 0.998, respectively. As both p-values are high, the null hypothesis that the observed and expected proportions are the same is accepted, i.e., the presented model (6) fits the given data well.
From the above results, it can be concluded that the model presented in (6) can provide a real-time prediction of the probability of the voltage sag duration based on the HF value of MHDU. Therefore, this hypothesis was confirmed, and the RTPP method could be applied.

F. DISCUSSION
The proposed RTPP method represents a binary classifier, where outputs are ''0,'' i.e., there is a high probability that sag will be ''momentary'' and ''1,'' i.e., there is a high probability that sag will be ''temporary''. However, the distribution of these solutions is unequal, and classification accuracy alone can be misleading. Therefore, the confusion matrix was chosen as an adequate technique for summarizing the performance of this classification algorithm [72].
The confusion matrix consists of rows representing the instances in a predicted class and columns representing the instances in the actual class. To fill the matrix, a probability of 0.5, for the trigger value, is chosen, i.e., if f (X = x) ≥ 0.5, then the sag with MHDU = x is predicted as temporary. The confusion matrix for the proposed model is presented in Table 3. It can be seen that there are only 2 false-positive cases (when momentary sag is predicted as temporary) and 2 false-negative cases (when temporary sag is predicted as momentary). Furthermore, three of those four false cases are in the uncertain zone, while one temporary is classified as momentary.
Considering all recorded signals, for only four (out of 95) duration probabilities could not be determined correctly (95.8% reliability). In Fig.1, they are marked with red hexagrams, with the sag type and number next to the mark.
Thus, it can be concluded that the reliability (success rate) of the algorithm across all recordings was 95.8%. Three out of four sags with incorrect estimations were in the uncertain zone, and for them, the model did not provide any prediction. These three sags were very short, with durations of 0.1 s, 0.05 s, and 0.05 s, respectively. The false estimation of the duration probability is for one temporary interruption because this fault contains multiple disturbances.

IV. THE RTTP ALGORITHM
The RTPP method is a real-time method; therefore it is suitable for implementation in the control system of power electronic converters connected to the grid (grid-tied inverter in the case of DER). A flow chart of the RTPP algorithm is presented in Fig. 3 as part of the control structure for dealing with voltage disturbances and LVRT requests.
The first step in all systems is to acquire three-phase voltage signals in digital form from the voltage sensors (transducers) in the grid. They are converted to a digital form by A/D conversion and sent to the power electronic converter inputs (Va, Vb, and Vc inputs in Fig. 3). The time required for preprocessing is normally not included in the voltage sag detection time, as it is usually much shorter. For example, the Hall voltage sensor response and A/D conversion time may be estimated to be less than 0.24 [ms].
The calculation procedure for the RTPP determination starts with the transformation of digital voltage signals from the time to the frequency domain using one of the well-known methods of time-frequency transform (DFT, RFFT, SFTF, or other). From the transformed signal (in the frequency domain), a specific set of low-order harmonic amplitudes are extracted (2 nd , 3 rd , 5 th , and 7 th harmonics), and HF (or HDU2357) is calculated. An HF pattern was formed (sampleby-sample) and its characteristics were determined.
The first monitored samples were treated using a neural network (NN) to obtain the early detection of voltage sags, as presented in [40]. If the sag is detected, then the disturbance mode is activated. The next step was to determine the maximum value of the HF or MHDU. The MHDU is typically reached in the first half-cycle of a fault (10 ms), therefore this parameter can be obtained in the 10-15 ms range. It can be used to predict two sag characteristics: magnitude and duration.
For magnitude (MoVS) prediction, the method presented in [46] was used. In this case, LVRT requirements can be fulfilled, and the DER may function in the grid supporting mode.
For the prediction of duration probability, which is the topic of this paper, the obtained value of MHDU=x is inserted in (6), and the probability of voltage sag duration is calculated (f (X = x)). The sag is then predicted as ''momentary'' or ''temporary'' using the S-curve shown in Fig.2.
The algorithm for RTPP determination and duration anticipation, as described above, is shown in Fig. 3. Together with other blocks, the method may be used as a pre-control block in an inverter control system for operation in cases of voltage sags.

V. APPLICATION OF THE RTTP METHOD A. APPLICATION OF RTPP FOR GRID-TIED INVERTER
The RTPP method, together with methods for voltage sag detection and MoVS prediction, are applied to the grid-tied  inverter control used for the DER (Fig. 4). The control system of the DER (including LVRT control) is modeled according to [73] and further modified with a pre-control block, which includes sags detection, MoVS prediction, and the RTPP subblocs, together with additional digital switches (as shown in Fig. 4).
In the detection block, the HF is calculated, and if sag is detected, the digital switch (SW1) is turned on. For voltage sag detection, smart algorithm based on HF using Recurrent neural network is used, as described in [40]. The maximum value of HDU2357 (MHDU) can be then determined. The value of MHDU is used in the RTPP algorithm; therefore, the voltage sag duration can be predicted to be momentary or temporary with corresponding probability. In addition, MoVS is predicted.
Based on the categorization results (momentary or temporary) and following the estimated MoVS values, further actions are taken regarding the operation in the LVRT mode. If the fault duration is predicted to be momentary and if MoVS is below the LVRT line, the disconnection will be delayed. Delay is set to be in the range of 0.3 s -1.5 s, depending on MoVS, and an extended LVRT control strategy can be determined, together with appropriate references for P and Q generation during the sag. In addition, information about the expected voltage levels during a fault is known after around 10 ms. Based on the requirements of IEEE Std. 1547-2018 and/or relevant grid-code, the minimum duration of LVRT support (before DG disconnection) can be approximately determined. The RTPP algorithm was designed to enable duration estimation at the beginning of the fault.

B. APPLICATION OF RTPP FOR GRID-TIED INVERTER IN THE MICROGRID HIL LABORATORY TESTBED
The complete method was further applied to the Microgrid Typhoon-HIL testbed in the Laboratory of Typhoon-HIL in Novi Sad, on the Typhoon HIL604 real-time emulator. The Typhoon HIL Control Center (THCC) program package was used together with the Schematic Editor and HIL SCADA subsystems. The power production during LVRT and the operation of the protective and switching equipment are controlled externally using an algorithm that includes the RTPP method.
The IEEE 33-bus grid was adapted with a 2 MW wind turbine. It consists of 33 bus bars (three-phase nodes) arranged on 4 main feeders. Each line connecting the nodes was 1 km long. The nodes are connected to 32 three-phase lines. The loads are connected to each of the 33 nodes. The source of the grid is a generator (component ''Grid'') that represents a balance node, a three-phase source, of the nominal voltage of 12.66 kV, which can generate the necessary amount of energy to sustain voltage. All other nodes (considering the power flow angle) are P/Q types. The total active power was 3.715 MW, whereas the total reactive power was 2.3 MVAr. The system is operating at a frequency of 60 Hz. A simplified diagram of the modified IEEE33-bus test grid is shown in Fig. 5.
A generic wind turbine, which is located within the microgrid subgroup in Library Explorer, was selected for the distributed power source in this system. The wind turbine is connected to the network between nodes 706 (node 6) and 726 (node 26) using a transformer of 0.4/12.66 kV/kV.
The Microgrid Typhoon-HIL testbed that is used for laboratory testing is shown in Fig. 6. The main part and the ''brain'' of the setup are Typhoon HIL 604 devices (three devices, marked in Fig. 6 with numbers 1, 4, and 6), which VOLUME 10, 2022 are used in the HIL emulation procedures. Below, there is the generator controller the Woodward easYgen-3500 (marked as 2). The SEL-751 relay (marked as 3), which is used to protect feeders in radial and tangled distribution networks is under it. Next in line, is another HIL 604 device (marked as 4), which is connected to the Typhoon HIL Connect device below it (marked as 5). The HIL Connect is used to interface the controller with the existing HIL system. Below is another HIL 604 (marked as 6), which is connected to a Typhoon HIL EPC Power controller (marked as 7).

C. APPLICATION OF RTPP IN THE FIELD OF PROTECTION
The method of voltage sags characterization using HF also has some advantages in the field of protection. The use of fast voltage sags detection in microprocessor relay operations has been reported [74].
The application of the proposed RTPP algorithm enables swift prediction in a period far below the normal reaction time of the protection relays. For example, the preprocessing stage (Hall voltage sensor response and A/D conversion) may be estimated to be less than 0.24 [ms], while the RTTP of additional 10-15 [ms]. However, false detection can lead to undesired disconnection of the DER or unselective/unwanted protection equipment operation. For example, the starting of induction motors and transformer energizing do not require the operation of a protection system, but they may be mistaken for sags caused by faults [46]. Therefore, it may have value in a modern microgrid for proper coordination between the protection system and LVRT requirements of the DERs [75]. In our simulation and experiment, we included the LVRT requirements, but relay tripping was not considered.

VI. TESTING OF THE PROPOSED METHOD
The proposed RTPP method is based on the HF feature of the voltage sag, i.e., on implementation of the sag's MHDU value. As mentioned previously, a set of real grid measurement results was used to support the method. Their HF parameters were determined and are presented in [46], where a strong correlation between the MHDU and MoVS was observed. In this paper, it is shown that by using the same measurement results and MHDU values, the voltage sag duration can be predicted with high probability. Complete tables with the test results are provided in the Appendix. In Tables 6 and 7 types of faults, the magnitude of voltage sag (MOVS), duration of fault, and maximum of harmonic footprint (MHDU) result are presented for momentary and sustained disturbances, respectively. In Table 8, the fits and statistical parameters for all real grid measurements are shown.

A. TEST SETUP -MODIFIED IEEE 13-BUS TEST GRID
To test the RTPP algorithm, real grid measurement data were applied to a distribution network model. The aim was to determine the benefits of the algorithm to the voltage supply of residential consumers. A modified IEEE13-bus test grid with a total of 3.25 MW of added DER generations (two PVs of 1.25 kW and one wind plant of 2 MW) and 0.464 kVA of residential area load was used. The IEEE13-bus test grid is a reduced model of an unbalanced real grid; therefore, with applied real-grid (in-field) measurement data, appropriate emulation of real grid cases may be achieved.
A simplified diagram of the modified IEEE13-bus test grid is shown in Fig. 7. A 2 MW wind plant consisting of four 0.5 MW wind turbines was added at a location 1 km from bus #633 (close to the distribution transformer). Wind turbines are doubly-fed induction generators (DFIGs) with an AC/DC/AC IGBT-based PWM voltage source converter (VSC). A 1 MW PV plant with a grid-tied inverter (DC/AC converter) was connected to the #652 bus through a 2 km line. A residential area with an added total of 0.25 MW roof-top PV systems and a total load capacity of 0.494 kVA was connected to the #634 bus.

B. TEST RESULTS WITH EMULATED REAL-GRID SAGS
The algorithm was tested using emulated real-grid measurement results, i.e., assuming that such voltage sags were recorded at the observed bus #634. A critical distance method was used to detect the location of faults in the IEEE 13-bus test grid [48], [49]. The faults are located around buses #632 (Fault Location 1, FL1) and #648 (FL2) and close to bus #692 (FL3), as indicated in Fig. 7. Therefore, the voltage sags at bus #634 resulting from faults in the grid have the same parameters (magnitude and duration) and time waveforms as the sags obtained from the recorded real-grid measurements (Fig. 1).
For testing, six voltage sags were selected corresponding to the ones shown in Fig. 1. They are outside or on the line of the German grid code minimum requirements for the LVRT, i.e., of such duration (0.41 s -1.55 s) and magnitude that they would probably cause a DER disconnection. The selected sags were regenerated using the detected fault locations in the test grid model and the appropriate fault parameters. For example, for fault #1, the parameters are fault  resistance 0.0001 , ground resistance 0.001 , fault/sag duration 0.41 s, and voltage sag MoVS 6.35%. The fault number, type, location, and parameters (duration and magnitude-MoVS (left number)) are presented in Table 4. Voltage sag #1 is shown in Fig. 8, sag #6 in Figs. 9 and 10, whereas the others are omitted because of lack of space. Simulation is building up at the start, and steady-state is reached at 0.05 s, so results from 0 s to 0.05 s may be discarded. The faults were generated at 0.09 s.
The values of the obtained MHDU and RTPP for all the six tested faults are presented in Table 4. The first (left) MoVS value is obtained when all DERs are disconnected (without RTPP), whereas the second (right) is obtained if they are active in the LVRT mode because of the RTPP usage. It can be observed that all MHDU values are below the limit of 37.7% (see Fig. 2); therefore, their RTPP is low, and all sags are anticipated to be momentary. Results show that for all six cases the RTPP algorithm successfully predicts duration as ''momentary'' with probability in the range from 85.47% for the fault (sag) #1, up to 99.63%, for the fault (sag) #3.
A graphical representation of the test results for fault #1 (SLG) is presented in Fig. 8. In Figs. 8 (a) and 8 (b), the line voltages at bus #634, with and without the integrated RTTP algorithm are presented, respectively. As shown in Fig. 8 (b), after 0.36 s the DERs are disconnected, following the LVRT minimum requirements from the German grid code. On the other hand, Fig. 8 (a) shows that the DER disconnections are prevented owing to the application of the proposed algorithm. In Fig. 8 (c), the harmonic footprint is presented with the important points indicated. The starting voltage sag's transient and ending transient can be distinguished. MHDU value (31.8%) can be observed, so the sag is predicted as ''momentary.'' Subsequently, the LVRT mode is activated, and the DERs support the grid; therefore, the voltage magnitudes during the sags at bus #634 are higher.
In Figs. 8 (d) and 8 (e), the active (P) and reactive (Q) power of the 2 MW wind farm with the proposed algorithm (blue  Table 4: (a) Voltage signals at a residential area in bus #634 without the RTPP algorithm; (b) Voltage signals at a residential area in bus #634 with the algorithm for the RTPP.  line) and without (red line) are presented, respectively. The wind farm was chosen because it has the highest power. The responses of PV plants were similar. The grid support criteria are defined by LVRT requirements, which may differ from country to country [22]. Usually, it is required to generate both P and Q if voltage during the sag is above 50% of the rated, and the Q only if it is less. It can be seen that the DER continues to support the grid, so voltages at bus #634 are better and DER operation is not disrupted.
In Fig. 9, another example of voltages in a residential area during a voltage sag in the grid is presented (single-phase to ground). The simulated sag, with a magnitude of 52.   other two phases, after DGs cut-out, an immediate increase in voltage can be observed, causing voltage swells in these phases.
In Fig. 10, the P and Q productions of the 2 MW wind power plant during a fault in the grid (fault #6 from Table 4 (d) that with the application of the algorithm presented in this paper, unnecessary tripping of the 2 MW wind power plant is avoided, and immediately after recovery of the voltage, normal power production is restored.
As MoVS during sag is 53%, according to [23], the production of active power during LVRT should be close to zero, and almost all delivered power should be reactive to support the voltage. The inverter control in the simulation was modeled following these recommendations.

C. TEST RESULTS WITH TYPHOON-HIL SETUP
The algorithm was further tested in a laboratory HIL environment using an advanced Typhoon-HIL Microgrid setup. For testing, five voltage sags were simulated in the modified IEEE 33-bus test grid. Faults were set to be located between buses #2 and #3 (Fig. 6). The measurements were performed at the interconnection point of the wind turbine to the grid. The results are presented in Table 5. It can be observed that all MHDU values are below the limit of 37.7% (see Fig. 2) and that all sags are anticipated to be momentary.
As shown in Fig. 11, the SCADA window is made up of a graphical representation of the IEEE 33 bus network, in which each node has a special sub-window showing the value of the phase voltage, as well as the phase shift relative to zero nodes. At the root of the network, there is a window that enables the setting of the grid component's voltage reference value and frequency as well as tracking these values, along with the active and reactive power of the network.
In Fig. 12, the line voltages for two cases, when the wind power plant control is with the integrated RTTP algorithm and without it, are presented. In Fig. 12 (a), it is shown that the DER disconnection is prevented by the application of the proposed algorithm. On the other hand, in Fig. 12 (b), it can be seen that after 1.25 s the DER is disconnected, following the LVRT minimum requirements from the German grid code, so another voltage transient followed by a further decrease in voltage can be observed. Figs. 12 (a) and 12 (b) are not aligned, and the starting time of recording of signals and the end time are not the same, as the testing is performed in real time. Similar to Fig. 12, the test results from the HIL testing are presented for Fault #5 in Fig. 13. In Fig. 13 (a) DER disconnection is prevented due to operation of the RTPP, while in Fig. 13 (b) it can be seen that after 0.72 s the DER is disconnected. Fig. 13 (c) shows the HDU2357 values during the voltage sag.
The results show that for four cases the RTPP algorithm successfully predicted the duration to be ''momentary'' with probability in the range from 90% for sag #5, up to 99.87% for sag #1, and sag #4 is successfully predicted to be ''temporary'' with the probability of 99%.

D. DISCUSSION
In Fig. 1, German grid code limits were added to indicate the minimum requirements for LVRT operation. They defined that a DER should support the grid within a certain period and then disconnect if a fault is not cleared. Considering the case of Fault #1 shown in Fig. 8, it should be at 0.25 s from the start of disturbance, that is, at 0.34 s (Fig. 8 (b)). Immediately after disconnection (0.36 s), the voltage magnitude dropped to zero ( Fig. 8 (b), blue line). In addition, after the DERs cut-off, an immediate increase in voltage can be observed, causing the voltage to swell in other non-faulted phases.
The voltage diagrams in Fig. 8 (a) show that with the implementation of RTPP, the sag duration is predicted as ''momentary''. This information is provided to the adapted control algorithm, which enables additional time for LVRT support -from 0.34 s to 0.5 s. At the 0.5 s, the voltage sag ends, voltage is recovered, and the DERs continue their normal operation without tripping.
Furthermore, similar results were obtained during laboratory testing with the Microgrid HIL testbed. The voltage     diagrams in Figs. 12 (a) and 12 (b) show similar behaviors to those presented in Figs. 8 (a) and 8 (b). After 1.25 s from the start of the disturbance (MoVS 71%), the DER is starting to disconnect from the grid (Fig. 12 (b)). Immediately after disconnection, the voltage magnitude decreased further to 43%. In Fig. 12 (a), it can be observed that with the RTTP algorithm the wind power plant continues the LVRT operation without tripping, and the voltage levels are stable throughout the entire sag.
Without this algorithm, a DER would be disconnected, and the grid support would be terminated. Therefore, some  amount of energy generation would be lost, reconnection would require additional time, and instability of the grid would be induced. Additionally, an improvement in LVRT support could be achieved through a faster response to the fault.
It can be seen that as a result of the RTPP application, the inverter control system can achieve a time gain for moving the operating point according to the LVRT requirements and starts supporting the grid. The time gain (TG) in [s] is the difference in the timeline between the instant when the duration of the fault is estimated with the RTPP and when the sag's end is recognized with a standard algorithm (RMS). This is presented in the last column of Table 4. Based on the five voltage sags presented in Table 4, the average time gain was 1.002 s.

VII. CONCLUSION
In this paper, a novel algorithm for real-time prediction of the probability of voltage sag duration based on HF is presented. Mathematical and statistical analysis of the maximum value (MHDU) of the HF relation to the probability of sag duration is provided with adequate verification using extensive real-grid measurement data. It was shown that a voltage sag duration can be predicted (with high probability) either as ''momentary'' (if MHDU<37.7%) or as ''temporary'' (if MHDU>41.2%). In between, is the uncertain zone (37.7%<MHDU<41.2%).
The algorithm is applied to the control of grid-tie converters used in DERs and tested with real grid measurement data in the IEEE 13-bus test grid by simulations, and in the IEEE 33-bus test grid using HIL Microgrid laboratory testbed. It is shown that unnecessary tripping of wind farms and PV plants can be avoided in cases of different types of momentary sags, with benefits to residential area customers. The benefits to customers are reflected in better voltages during the faults and in avoiding unnecessary instability.
This research was based on a large set of real grid recordings and their analyses. Two different voltage levels guaranteed the diversity of the data. However, a specific S-curve that shows high reliability may differ. As this is a onedimensional problem, it is possible to improve the reliability using neural networks to learn and adjust the probabilistic model for specific applications and under real operating conditions. APPENDIX Tables 6 and 7 present a detailed description of the faults obtained by measurements in a real grid, which were used in the research. Table 8 shows the fits and diagnostics for all observations.