1. Introduction

The increase in the use of renewable energy sources combined with energy efficiency measures has reduced consumption and imports of fossil fuels and has consequently contributed to increase the representativeness of renewable energies in final energy consumption. The overall capacity of all wind turbines installed worldwide by the end of 2017 reached 539 GW ^[1]. The annual onshore and offshore wind installations achieved, respectively, 12.5 GW and 3.1 GW in EU. There are differences between onshore and offshore exploitation of wind energy, mostly due to the: the wind speed is relatively more stable at offshore than at onshore; the wind speed at offshore is about 20% higher than at onshore. Both these differences imply that the utilization rate is higher at offshore than at onshore. Economically wind energy conversion systems at offshore tend to have larger capacity and benefit from lager areas than the limited and mostly occupied areas at onshore. But, some issues must be regarded with caution, for instance, regarding preservation of fishing, boat navigation, ocean habitat and species, marine conservation zone, landscape view from land to sea and vice versa. Higher reliability, resistant to corrosion due to the difficulties and costs of maintenance or repairing at offshore and consideration of action of sea waves on wind energy conversion systems are needed. Turbines need to be designed to take advantage of the higher level of the wind speed, i.e., the design must consider a higher value for the cut-out wind speed: the wind speed at which the wind turbine shuts down to prevent mechanical damage of the structure. Foundations work at offshore are more expensive, than on onshore where the process uses either rods drilled into the ground or a reinforced concrete pad set in the ground. Turbines need to be designed to take advantage of the higher level of the wind speed at offshore, i.e., the design must consider a higher value for the cut-out wind speed: The wind speed at which the wind turbine shuts down to prevent mechanical damage of the structure. Foundations work at offshore are more expensive, than on onshore where the process uses either rods drilled into the ground or a reinforced concrete pad set in the ground.

The technology applied on wind energy conversion systems having variable-speed and variable-pitch (VSVP) has been reported in literature due to advantageous energy capturing ^[2,3]. Also, amongst the electrical generators available for equipping wind turbines having VSVP, doubly-fed induction generator (DFIG) stands out in industry and became the mainstream choice ^[4,5].

Control system architectures are necessary to improve the quality of electrical energy delivered into the electric grid. Pitch control is widely used to ensure the best performance during the capture of energy under all operational wind scenarios ^[6,7,8]. Control strategies have to deal with unforeseen variables such as wind speed variability and intermittence, to achieve the goal of an overall acceptable performance. Various control techniques are used by researchers in the architecture of control strategies, amongst a wide variety of control types, some are mentioned here: Robust multivariable gain-scheduled control ^[9], fuzzy proportional integral ^[10], adaptive linear quadratic control ^[11], data-driven and adaptive PI control ^[12], sliding-mode control ^[13].

A suitable design for control and supervision under the unpredictable character of the wind energy is challenge in order to maximize power output. Thus, the controller needs to have access to the variables generator speed and pitch angle. The performance during the capture of wind energy is a prominent issue and has been dealt not only by research on turbine design and advanced control strategy, but also by condition monitoring using a supervisor. Some works addressed this issue, such as: Classification of wind turbine condition monitoring methods and techniques with a focus on trends and future challenges ^[14], supervisory control and data acquisition system ^[15]. The supervisory control theory is behind the architecture proposed in this paper and is suitable for control application with event-based operations as can be seen in previous works ^[16,17].

This paper presents the following contributions: A two-level control architecture, consisting in a supervisory level (event-based supervisor) and an execution level (LQG and MPC controllers); a study and implementation of adaptive controllers and predictive controllers; a study of an event-based supervisor, having the objective of determining the most suitable operational state among the possible states; a simulation of the system with the distinct types of controllers with and without the presence of an event-based supervisor; and an evaluation of the system performance driven from the assessment of the results obtained with the action of each controller. In the overview of this paper, Section 2 presents the onshore energy conversion system modelling, Section 3 presents the control and supervision strategy applied to the system, namely the LQG, MPC in the execution level, an event-based supervisor in the supervision level and is also defined the operating regions according to the wind speed. Section 4 presents the simulation results as well the performance of the system without and with supervisor. Conclusions remarks are given in Section 5.

---

2. Onshore energy conversion system modelling

A mathematical model representing an appropriate and thorough dynamic of an onshore energy conversion system is developed in this section. The development of the mathematical model is based on the benchmark model developed in ^[19], whose electric power output is 4.8 MW. The proposed control and supervision hierarchical structure is represented by a block diagram model composed by the following functional systems: Mechanical system (blade-pitch angle and two-mass drive train), controllers (LQG and MPC), electrical system (generator and two-level converter) and event-based supervisor. The block diagram model shown in Figure 1, is composed by the following variables: v_w is the wind speed in m/s, τ_g and τ_r are the generator and turbine rotor torques in Nm, ω_g and ω_r are the generator and turbine rotor speed in rad/s, β_r is the pitch angle in degrees, P_g and P_r are the generator and reference power in MW. The r and m subscripts designate respectively reference and measurements values.

---

2.1. Blade and pitch modelling

This model combines the aerodynamic with blade and pitch models. The aerodynamic torque is given by:

${\tau _r}(t) = \frac{{\rho \pi {R^3}{C_p}\left( {\lambda (t), \beta (t)} \right){v_w}{{(t)}^2}}}{2}$

(1)

Where R is the radius of the blades, $\rho$ is the air density, ${C_p}\left( {\lambda (t), \beta (t)} \right)$ is the power coefficient, which is a function of the tip speed ratio $\lambda (t)$ and the pitch angle $\beta (t)$. The tip speed ratio is given by:

$\lambda (t) = \frac{{{\omega _r}(t)R}}{{{v_w}(t)}}$

(2)

From Eq 2, variations of the wind speed can lead to two consequences:

a) If the turbine rotor speed is held constant, then λ(t) will change, leading to a consequent change in the power coefficient;

b) If the turbine rotor speed is suitably adjusted, then λ(t) can be held at a reference point and as a result power coefficient can be kept at a desired value.

The power coefficient of a wind turbine using pitch control ^{[20,21,22,23]} is given by:

$\left\{ \begin{array}{l} {C_p}(\lambda , \;\beta ) = 0.73\, {\kern 1pt} ({\kern 1pt} {\kern 1pt} \frac{{151}}{{{\lambda _i}}} - 0.58\, \beta - 0.002\beta {\kern 1pt} {{\kern 1pt} ^{2.14}} - 13.2{\kern 1pt} {\kern 1pt} {\kern 1pt} ){\kern 1pt} {\kern 1pt} {\kern 1pt} \, {e^{\frac{{ - 18.4}}{{{\lambda _i}}}}}\\ {\lambda _i} = \frac{1}{{\frac{1}{{(\lambda - 0.02\, \beta )}} - \frac{{0.003}}{{({\beta ^3} + 1)}}}} \end{array} \right.$

(3)

The pitch angle can be modelled as a second order differential equation given by:

$\ddot \beta (t) = - 2\xi {\omega _n}(t)\dot \beta (t) - \omega _n^2\beta (t) + \omega _n^2{\beta _r}(t)$

(4)

---

2.2. Drive-train modelling

The mechanical drive train is composed by two shafts, a low-speed shaft from the rotor side, a high-speed shaft from the generator side, connected by a transmission equipped with a gear box. The two-mass drive train is shown in Figure 2.

The mathematical model to describe the dynamic of the drive train is a two-mass model. On the rotor side, the mass J_r concentrates the inertia of the hub, turbine blades and low-speed shaft inertia. On the generator side, J_g concentrates the generator inertia and high-speed shaft. The low-speed shaft and high-speed shaft are connected by a gear box with a ratio N_g. The torsion shaft stiffness is K_dt and torsion shaft damping is B_dt. The angular deviation due to the damping and stiffness coefficients between the turbine and the generator is ${\theta _\Delta }(t)$. The linear model for the drive train model is given by:

${J_r}{\dot \omega _r}(t) = {\tau _r}(t) + \frac{{{B_{dt}}}}{{{N_g}}}{\omega _g}(t) - {K_{dt}}{\theta _\Delta }(t) - \left( {{B_{dt}} + {B_r}} \right){\omega _r}(t)$

(5)

${J_g}{\dot \omega _g}(t) = \frac{{{K_{dt}}}}{{{N_g}}}{\theta _\Delta }(t) + \frac{{{B_{dt}}}}{{{N_g}}}{\omega _r}(t) - (\frac{{{B_{dt}}}}{{N_g^2}} + {B_g}){\omega _g}(t) - {\tau _g}(t)$

(6)

${\dot \theta _\Delta }(t) = {\omega _r}(t) - \frac{1}{{{N_g}}}{\omega _g}(t)$

(7)

---

2.3. Generator and converter modelling

The model for the generator and the power converter is described by the state equation given by:

${\dot \tau _g}(t) = - {\alpha _{gc}}{\tau _g}(t) + {\alpha _{gc}}{\tau _{g, r}}(t)$

(8)

Where ${\alpha _{gc}}$ is a first order constant and ${\tau _{g, r}}$ is the reference torque of the generator. The power of electric generator is given by:

${P_g}(t) = {\eta _g}{\omega _g}(t){\tau _g}(t)$

(9)

Where ${\eta _g}$ is the efficiency of the generator.

---

3. Control and supervision strategy

The onshore energy conversion system components, the wind speed and the wind turbulence, the tip speed ratio and the pitch angle are some of the variables considered in the control design to achieve the rated power with an acceptable overall performance. In this paper two distinct controller are considered, adaptive and predictive. The control and supervision strategy are based on the switching from Region Ⅱ (constant-pitch, generator torque control) to Region Ⅲ (variable-pitch control, constant torque). Region Ⅱ (wind speed is between 5 m/s and 13 m/s) and Region Ⅲ (wind speed is between 14 m/s and 22.5 m/s) can be seen at the regions of power operation as shown in Figure 3.

Figure 3 shows that in Region Ⅱ, the control objective is to extract optimal power reducing oscillation on the drive train, thus torque control is applied. The pitch angle and tip speed ratio are set to their optimal value to achieve this control objective and the rotor speed must be adjusted to the reference rotor speed given by:

${\omega _{ref}}(t) = \frac{{{\lambda _{opt}}(t){v_w}(t)}}{R}$

(10)

The optimal electric generator torque is given by:

$\left\{ \begin{array}{l} {\tau _{g, r}}(k) = {K_{opt}}{\left( {\frac{{{\omega _{ref}}(k)}}{{{N_g}}}} \right)^2}\\ {K_{opt}} = \frac{1}{2}\rho A{R^3}\frac{{C{p_{\max }}}}{{\lambda _{opt}^3}} \end{array} \right.$

(11)

Where ${\lambda _{opt}}$ is the optimal point for achieving maximum C_P and A is the area swept by the blades.

Figure 3 shows that in Region Ⅲ, the control objective is to control the output generator power through pitch control. To achieve this control objective, the pitch angle varies according the error between nominal speed and the output generator speed. The pitch reference and generator torque reference should be adjusted at the same time and the latter is given by:

${\tau _{g, r}}(k) = \frac{{{P_r}(k)}}{{{\eta _g}{\omega _g}(k)}}$

(12)

Switching between Region Ⅱ and Region Ⅲ is determined by the electric generator power, P_g, and the generator speed, ω_g, the control mode switching is performed according to the following conditions: Switching from Region Ⅱ to Region Ⅲ is enabled if ${P_g} > {P_r}\, \, \, and\, \, \, {\omega _g} > {\omega _{nom}}$ and switching back from Region Ⅲ to Region Ⅱ is enable if ${\omega _g} < {\omega _{nom}} - {\omega _\Delta }$. The non-null threshold ${\omega _\Delta }$ is required to enable preventing sudden switches from Region Ⅲ and Region Ⅱ.

---

3.1. Adaptive linear quadratic gaussian controller

Discrete adaptive linear quadratic gaussian control (LQG) deals with unpredicted variables and the adaptation depends on estimation of parameter $\hat \theta \left( k \right)$. Using recursive least squares algorithm, polynomials parameters from A₁(z⁻¹) and B₁(z⁻¹) can be estimated, obtaining parameters $\hat \theta \left( k \right) = \left[ {{{\hat a}_{11}}\;{{\hat a}_{12}}\;{{\hat b}_{11}}\;{{\hat b}_{12}}} \right]$. The dynamic development of this controller is presented in ^[24]. The dynamical system of the onshore wind turbine benchmark around a nominal set-point r(k) is represented by an ARX model, thus the linear quadratic control transfer function is given by:

$\frac{{Y(z)}}{{U(z)}} = \frac{{\left( {{b_{11}}{z^{ - 1}} + {b_{12}}{z^{ - 2}}} \right)}}{{\left( {1 + {a_{11}}{z^{ - 1}} + {a_{12}}{z^{ - 2}}} \right)}}$

(13)

Where the Z-transform of the system output and control input are respectively Y(z) and U(z). The performance index H is given by:

$H\left( k \right) = {\left( {Py(k + d) - Qr(k)} \right)^2} + {\left( {\chi u(k)} \right)^2}$

(14)

Where P and Q assume a unit value, the scalar $\chi = 0.4$ is chosen to achieve an acceptable closed-loop performance and optimal control u(k) minimizes performance index.

The LQG pitch control structure is shown in Figure 4.

The polynomials F(z⁻¹) and G(z⁻¹) can be determined by solving the Diophantine equation. The optimal control that minimizes the performance index is given by:

$u(k) = \frac{{\mathop {{b_{11}}}\limits^ \wedge }}{{{{\mathop {{b_{11}}}\limits^ \wedge }^2} + {\chi ^2}}}\left[ {\left( {\frac{{{\rho ^2} - \mathop {{b_{12}}}\limits^ \wedge \mathop {{b_{11}}}\limits^ \wedge }}{{\mathop {{b_{11}}}\limits^ \wedge }}} \right)u(k - 1) + \mathop {{a_{11}}}\limits^ \wedge y(k) + \mathop {{a_{12}}}\limits^ \wedge y(k - 1) + r(k)} \right]$

(15)

---

3.2. Model predictive controller

Model predictive controller ^[25] use the states of the blades and pitch model, drive train and nonlinear aerodynamics as input. By optimizing the best output according to the operating region and control objective, the model predictive controller generates reference values for β(k) and τ_{g, r} (k). From Eq 1 to Eq 9 a linear state space model representing the onshore energy conversion system dynamics at certain operating wind speed is given by:

$\left\{ \begin{array}{l} {\mathit{\boldsymbol{x}}}\left( {k + 1} \right) = {\mathit{\boldsymbol{Ax}}}\left( k \right) + {\mathit{\boldsymbol{Bu}}}\left( k \right)\\ {\mathit{\boldsymbol{y}}}\left( k \right) = {\mathit{\boldsymbol{Cx}}}\left( k \right) \end{array} \right.$

(16)

Where $\mathit{\boldsymbol{x}} = {\left[ {{\omega _r}\;{\omega _g}\;{\theta _\Delta }\;{\tau _g}\;\beta } \right]^T} \in {\mathbb{R}^5}$ is the state vector, $\mathit{\boldsymbol{u}} = {\left[ {{\tau _{g, r}}{\beta _r}} \right]^T} \in {\mathbb{R}^2}$ is the control input and $\mathit{\boldsymbol{y}} = {\left[ {{\omega _g}\;{{\text{P}}_g}} \right]^T} \in {\mathbb{R}^2}$ is the measured output. One important purpose in this region is to maintain the generator speed very close to the nominal speed with the least stress on the generator and pitch control system. To accomplish this objective, one must follow the mathematical formulation for the objective function, i.e., the cost function is given by:

$\mathop {\min }_{u\left( k \right) \cdots u\left( {k + Np - 1} \right)} J\left( k \right) = \underbrace {\mathop \sum _{j = 1}^{Np} {{\left[ {r\left( {k + j} \right) - \hat y\left( {k + j|k} \right)} \right]}^T}Q\left( j \right)\left[ {r\left( {k + j} \right) - \hat y\left( {k + j|k} \right)} \right]}_{Quadratic\;Error} + \underbrace {\mathop \sum _{j = 0}^{Np - 1} {{\left[ {u\left( {k + j|k} \right)} \right]}^T}R\left( j \right)\left[ {u\left( {k + j|k} \right)} \right]}_{Control\;Effort}$

(17)

In Eq 17 the cost function is associated with the prediction horizon, given by the periods j with j = 1, …, Np. The first term is a quadratic functional, weighted by matrix Q, of the difference between reference $r(k + j)$ and predicted output $\hat y(\left. {k + j} \right|k)$ and is used for penalization of the error. Also, the second term is a quadratic functional, weighted by matrix R, measuring the control effort $u(\left. {k + j} \right|k)$.

The MPC pitch control structure is shown in Figure 5.

The optimization problem consists in the minimization of the cost function subjected to the following constraints: ${\beta _{\min }} \leqslant {\beta _r} \leqslant {\beta _{\max }}$, $\left( {{\omega _{nom}} - {\omega _\Delta }} \right) \leqslant {\omega _g} \leqslant \left( {{\omega _{nom}} + {\omega _\Delta }} \right)$ and also ${\nu _{\min }} \leqslant {\nu _w} \leqslant {\nu _{\max }}$. Where ω_nom is the nominal generator speed and ${\omega _\Delta }$ is a small offset that introduce some hysteresis in the switching scheme avoiding frequent switching from control modes.

---

3.3. Supervisor

The onshore energy conversion system event-based supervisor is based on finite state machine (FSM) also known by finite-state automaton or state machine. This machine is a mathematical model of computation to conveniently schedule the operational states of a process. The states are such as: start, production and stop. The machine is in only one state at a time and can change from one state to another when initiated by a triggering event or condition called a transition. The machine is a set of: behavior states, transitions that define the transit between states and rules or conditions to be fulfilled in order to enable a transition. The specification of a machine is given by the list of states and the triggering condition for transitions. The FSM is described by a quintuple model ^[26] given by:

$FSM = (\Sigma , {Q_f}, \partial , {q_0}, {Q_m})$

(18)

In Eq 18 the Σ is the input alphabet, ${Q_f}$ is the finite set of states, $\partial = \sum \times \, Q$ is the transition function, ${q_0}$ is the initial state and ${Q_m}$ is the set of final states. The event-based supervisor deterministic version of an FSM used in this paper has the state transition diagram shown by Figure 6.

In Figure 6, the operational states are park, start, generating and brake, typifying the regions of power operation as follow:

a) In the Park state, the wind turbine is in shutdown mode and the generator is disconnected from the electric grid. This state corresponds to Region I.

b) In the Start state, the wind speed must be above the cut-in speed (v_min). The wind turbine should be started to capture energy from the wind power. The generator is connected to the electric grid, but not necessarily at rated power in the majority of the operation in Region Ⅱ. From this state and depending on the accuracy of the transitions (wind speed and generator speed), one can enter into the generating state or into the brake state. This state corresponds to Region Ⅱ.

c) In the Generating state, also known as power production, the wind speed must be within the rated wind speed (v_rated) and the cut-out speed (v_max). The generator is connected to the electric grid at rated power, by conveniently curtailment of the capture of kinetic energy from the wind. This state corresponds to Region Ⅲ.

d) In the Brake state, also known as slow down state, the pitch angle is adjusted in order to decrease the rotation of the wind turbine. From this state and according to the value of the conditions, one can enter into the start state or into the park state.

---

4. Results

The mathematical model of the WECS equipped with a two-level power converter topology is used in the simulations. The simulations in order to compare the performance without or with the supervisor are carried out in Matlab/Simulink. The time horizon considered in the simulations is of 4500 s, and the sampling time Ts is of 0.01 s. The WECS has a rated electric power of 4.8 MW for a wind speed ${v_w}$ of 13 m/s, where ω_nom is 162 rad/s, ${\omega _\Delta }$ is 15 rad/s and ${\omega _{\max }}$ is 200 rad/s ^{[19,27,28,29]}. The series for the wind speed used in the simulations are the ones given in ^[19], having amplitudes in the range of 7.5 m/s to 22.5 m/s between Region Ⅱ and Region Ⅲ, but an addition of white noise to the wind speed is considered in order to have an augment effect of intermittence of wind energy. The wind speed with the addition of white noise is shown in Figure 7.

---

4.1. Without supervisor

The following figures represent the simulations of the controllers, LQG and MPC, without the presence of the event-based supervisor. The electric output and reference power of the LQG and the MPC controllers are shown in Figures 8a, b.

In Figure 8 both controllers, the electric output power follows the reference power with some oscillations around the reference power having some decreases in the output power. Major strengths of MPC are abilities to handle multivariable interactions and operating constraints in systematic manner. MPC is formulated as a constrained optimization problem, which is solved on-line repeatedly by carrying out model-based forecasting over a moving window of time. More importantly, MPC facilitates optimal control systems with unequal number of manipulated inputs and measured outputs. LQG provides a systematic approach to designing a control law for linear multi-variable systems. However, some difficulties associated with the classical LQG formulation is inability to handle operating constraints explicitly. Operating constraints, such as limits on manipulated inputs or on their rates of change, limits on controlled outputs arising out of product quality or safety considerations, are commonly encountered in any control application.

MPC can be viewed as modified versions of LQ (or LQG) formulation, which can deal with the operating constraints in a systematic manner. Origins of MPC can be traced to the classical linear quadratic optimal control (LQOC) formulation. These techniques, although have different mechanisms, are similar, thus, giving similar outputs.

The LQG controller and the MPC controller pitch angles are shown in Figures 9a, b.

In Figure 9 LQG controller provides pitch angle variations between 7 degrees and 18 degrees. MPC controller provides pitch angle variations between 10 degrees and 65 degrees.

The switching between Region Ⅱ and Region Ⅲ is shown in Figures 10a, b.

In Figure 10 Region Ⅱ (0) corresponds to startup operational mode and Region Ⅲ (1) corresponds to generation operational mode. Both controllers present, with some frequency, switching between Region Ⅱ and Region Ⅲ.

---

4.2. With supervisor

The following figures represent the simulations of the controllers, LQG and MPC, with the presence of the event-based supervisor. The electric output and reference power of LQG and the MPC controllers are shown in Figures 11a, b.

In Figure 11 MPC controller allows a smoother response around reference power, when compared to the one obtained with the LQG controller, however presents several peaks.

The LQG controller and the MPC controller pitch angles are shown in Figures 12a, b.

In Figure 12 LQG controller pitch angle variation only occurs in a small interval of time reaching the maximum of 40 degrees. LQG controller contributes to less variation on the pitch angle. MPC controller provides pitch angle variations between 15 degrees and 25 degrees having one peak above 30 degrees.

The switching between Region Ⅱ and Region Ⅲ is shown in Figures 13a, b.

In Figure 13 the transition between the two operational states is due to the need of sustaining the electric output at the rated power. The LQG controller allows a less frequent switching between Region Ⅱ and Region Ⅲ, with permanence in Region Ⅲ while MPC controller presents few switching between Region Ⅱ and Region Ⅲ.

---

4.3. Performance assessment

The metrics used in the evaluation of the performance of the controller are the integral of time multiplied by the absolute value of the error (ITAE) and the integral of the square value (ISV). The ITAE is given by:

$ITAE = \int_0^{{t_f}} {t\left| {e(t)} \right|} dt$

(19)

And the ISV is given by:

$ISV = \int_0^{{t_f}} {{u^2}(t)} dt$

(20)

The control performances for the WECS without or with supervisor control are summarized in Table 1.

In Table 1 considering the values obtained without the event-based supervisor, the error between the electrical output and the reference power is smaller with the LQG controller, although with the MPC the error is close to that presented by the LQG controller. As far as the control effort is concerned, the LQG controller has a smaller error when compared to the MPC controller, consuming less energy. Considering the values obtained under the event-based supervisor, both controllers accomplish equivalent performance with respect to the error between electric output and reference power. Regarding the control action effort, it is observed that the lowest effort is obtained with the LQG controller while the MPC requires more energy to maintain the output power around reference power.

---

5. Conclusion

An onshore energy conversion system connected to the electric grid and under an event-based supervisor control is presented. The control and supervision strategy proposed for the onshore energy conversion system is based on a hierarchical structure with two levels. LQG and MPC controllers are in the execution level whereas the event-based supervisor is placed in the supervision level.

The closed loop response between LQG and MPC without supervisor is quite similar, with abrupt variations around reference power. The pitch angle variations achieve lower amplitude values with the LQG controller. In terms of performance, LQG controller outperforms the latter at expense of lower control effort.

The closed loop response with supervisor remains similar as the above one, with smoother variations around reference power. The pitch angle variations with LQG only occurs in a small interval of time reaching the maximum of 40 degrees. In terms of performance, LQG controller maintains superiority over MPC.

---

Acknowledgments

This work is funded by: European Union through the European Regional Development Fund, included in the COMPETE 2020 (Operational Program Competitiveness and Internationalization) through the ICT project (UID/GEO/04683/2013) with the reference POCI010145FEDER007690; Portuguese Funds through the Foundation for Science and Technology-FCT under the project LAETA, reference UID/EMS/50022/2013; Portuguese Foundation for Science and Technology (FCT) under Project UID/EEA/04131/2013.

---

Conflict of interest

All authors declare no conflicts of interest in this paper.

---