An Approach to Determine the Weibull Parameters for Wind Energy Analysis: the Case of Galicia (spain)

The Weibull probability density function (PDF) has mostly been used to fit wind speed distributions for wind energy applications. The goodness of fit of the results depends on the estimation method that was used and the wind type of the analyzed area. In this paper, a study on a particular area (Galicia) was performed to test the performance of several fitting methods. The goodness of fit was evaluated by well-known indicators that use the wind speed or the available wind power density. However, energy production must be a critical parameter in wind energy applications. Hence, a fitting method that accounts for the power density distribution is proposed. To highlight the usefulness of this method, indicators that use energy production values are also presented.


Introduction
In recent years, the use of wind energy has been continuously growing, even at double-digit rates in several countries.Spain has the fourth largest wind capacity in the world rankings, with 21,674 MW installed in 2012 [1].In northwest Spain, Galicia represents approximately 15% of this capacity with 3311 MW, and there are plans in development to install more than 3300 MW.Furthermore, the Galician installed power density was 11.2 MW/100 km 2 in 2012, which is greater than the mean value of Spain (4.3 MW/100 km 2 ) and also greater than that in Denmark (9.2 MW/100 km 2 ), Germany (8.8 MW/100 km 2 ) and The Netherlands (5.5 MW/100 km 2 ).

OPEN ACCESS
In this context, characterizing the wind speed at a specific location or area is extremely important.This task is complex due to the random nature of the wind, which does not exactly follow any known statistical distribution [2]; this behavior is reflected, to a certain extent, in the power delivered by a wind turbine generator (WTG) [3][4][5].Therefore, the best way to characterize the wind speed at a specific location is to perform "in-situ" measurements, which should last several years.However, in certain studies, the use of known probability density functions (PDFs) [6][7][8][9] is useful or unavoidable, for example, in commercial software packages that specialize in wind energy resources [10], national wind energy resource atlases [11,12], international standards [13], simulations of WTG behaviors [14,15], development of site-matching approaches [16,17], etc.
Difficulty arises in choosing the best PDF that fits the wind speed distribution [6,8].Although there are several PDFs for this purpose [18][19][20], the PDFs most commonly used by researchers that study wind characteristic at wind sites are the Weibull and Rayleigh functions, which appear to be related to the nature of the wind in certain conditions [2,21].
Weibull parameters are typically obtained using well-known estimation methods, e.g., maximum likelihood, and the goodness of the resulting fits are evaluated by several indicators, e.g., R 2 [6,[18][19][20]22].In wind energy applications, evaluating the energy that can be produced in a certain area is one of the most important results associated with the estimation process [23].However, typical PDF fitting methods and fitness indicators do not specifically consider the way energy is produced by WTGs, i.e., their power curves.For this reason, a fitting method that takes into account the typical behavior of WTGs is proposed.Additionally, a set of fitness indicators that considers the performance of a WTG is also presented.
To evaluate the proposed fitting method and the proposed fitness indicators, data from several weather stations distributed around Galicia (northwest Spain) were used in this study.As background work related with this paper, certain studies that have had similar objectives must be emphasized.Carta et al. [6] presented an extensive review of typical PDF distributions and parameter estimation methods using wind data from the Canary Islands.Celik [18] and Akdağ et al. [24] specifically studied the estimation of the Weibull parameters in an area using different methods, where the results were evaluated using fitness indicators on the wind speed distribution and errors in the available wind energy.Chang [19], Seguro et al. [20], Stevens et al. [25], and Carta et al. [26] analyzed the relationship between WTG energy production and the fitness of a Weibull PDF.In these studies, WTG power curves were modeled by different equations using a reduced set of commercial machines as a reference.In this case, the error in the energy production was included as a fitness indicator but was not been included in the estimation methods of the PDF parameters.

Wind Speed Data in Galicia (Spain)
The analyzed region in this paper is Galicia, a region in northwest Spain, which is located in one of the windiest areas of Europe [27].The Galician Meteorological Service, Meteogalicia, which is an entity that depends on the local government, publishes data from all the weather stations distributed around Galicia on its web page [28].The data from these stations were used in this paper; however, only the windiest sites were selected to be studied.Wind speed data from the weather station in the Sotavento Experimental Wind Park [29,30], with measurements taken at heights of 20 m and 40 m, are also included in this paper.The 29 wind sites that were selected are listed in Table 1 and shown in Figure 1.

Wind Speed Distributions
Wind speed is typically measured at weather stations at a height of 10 m every 10 m.The resulting wind speed series can be represented as: [ ] where n is the number of data points.
The basic representation of wind speed data is the histogram.The most common form of the histogram is obtained by splitting the range of data into equally sized bins, called classes.Each class is represented by the middle value of the bin.Therefore, each bin b j with ∆v width has associated a relative frequency: where n j is the number of data points that falls inside the bin represented by the wind speed; v j and fr j is the relative frequency associated with bin "j".With this definition, the following relationships are fulfilled: n n; fr where N is the number of bins.Finally, the chosen bin width ∆v for the histogram can affect the fitting results [19,26].However, the typical values used for wind energy analysis are 0.5 m/s or 1 m/s [13,19,26].In this paper, following the international standards, the value of ∆v = 0.5 m/s was selected.
The most widely used PDF to fit wind data is the Weibull distribution, which is defined as [2,21]: ( ) where k is the unitless shape parameter, and c is the scale parameter in m/s.When k is equal to 2, another commonly adopted PDF is obtained, the Rayleigh distribution; this distribution is used in the studies in the international standard IEC 61400-12-1 [13].

Performance of Wind Energy Conversion Systems
The available power of the wind that crosses the rotor of a WTG is [9,31]: ( ) where p w (v) is the power associated with wind speed v; A is the rotor area; and ρ is the air density (typically 1225 kg/m 3 [13]).This power relates to the power generated by a WTG by means of the power coefficient: where p(v) is the power generated by the WTG; and c p is the power coefficient, which depends on the blade design, tip angle and the relationship between the rotor speed and wind speed.Its maximum theoretical value, known as the Betz limit, is 16/27 (≈0.593).However, this value is not achievable with real WTGs; typically, its maximum value is approximately 0.5 [32].The power coefficient can be obtained from manufacturer data where, apart from the aerodynamic behavior of the blades, the mechanical and electrical losses are considered.When evaluating the available energy at a wind site, the following function is used: This function is called the wind power density distribution and represents the distribution of wind energy at different wind speeds per unit of time and rotor area (W/m 2 ).For a specific wind site, it can be obtained from Equation ( 5): Therefore, the total wind power density E w is: ( ) ( ) When a Weibull PDF is considered, the following equation can be used [31]: ( ) where Γ is the Gamma function.The influence of the k and c parameters on the E w values will be discussed in Section 5.2.

Power Curves
The energy production of a WTG can be obtained by means of its power curve, where the relationship between the wind speed and the delivered power is established, and can be expressed by the following (see Figure 2) [32]: where p(v) is the electric power; v ci is the cut-in wind speed; v co is the cut-out wind speed; v r is the rated wind speed; P r is the rated power; and q(v) is a non-linear relationship between the power and wind speed.Several expressions can be used to represent the non-linear part of the power curve q(v) [12,33,34]; however, for the sake of simplicity, the cubic equation is typically used [32,[35][36][37]: where C p,eq is a constant equivalent to the power coefficient.Using the equation above, the relationship between the rated power and rated wind speed is 3 1 r , e q 2 ρ p r P AC v = .

Energy Evaluation
To evaluate the energy production, the distribution of the energy generated by a WTG at different wind speeds per unit of time and rotor area is considered, which is called the power density distribution e(v) in W/m 2 and is defined as: Taking into account Equation ( 11), the following equation can be used: The total power density E at a specific wind site and WTG can be obtained from: As a particular case, the power density E can be obtained analytically when the wind is represented by a Weibull PDF, and the WTG is modeled using a cubic power curve, see Equations ( 11) and (12).In this case, the following equation can be used [36,37]: where Γ is the gamma function, and γ is the incomplete gamma function.This equation is more complex than that used for the wind power density E w , shown in Equation (10).Apart from the inclusion of the terms ρ and C p,eq , this complexity is primarily related with the power curve behavior at the rated wind speed v r .

Weibull Parameters and Energy
The impact of Weibull parameters, k and c, on the energy production can be analyzed using Equation (16).With this purpose, a set of wind power density E w and power density E curves are shown in Figure 3, which were calculated using different scale and shape parameters (c and k).As shown in this figure, the relationship between the wind power density E w with respect to the Weibull parameters, c and k, is clear: higher E w values were obtained as c increased or k decreased.However, this statement is untrue when power density E is analyzed because for certain c and k values, the trend behavior is inverted.In conclusion, the power density uncertainty cannot be associated with any variation interval of k and c values.Therefore, the fitting methods must be analyzed separately due to the error related to their use, as shown in the following sections.

Power Curve Parameters
To characterize wind power curves, a database with WTG parameters [32] was used to obtain the relative frequencies and therefore, also obtain the distribution of the cut-in, cut-out and rated wind speeds, as shown in Figure 4.These values were used to define the cubic part of the power curve.After this analysis, the typical values for the cut-in wind speed (2.5-4.5 m/s), cut-out wind (20-35 m/s, typically 25 m/s) and rated wind speed (10-17 m/s) were obtained.The rated wind speed has a greater effect than the cut-in wind speed or cut-out wind speed from the point of view of energy production [25,34,38], which can be verified with Figure 5, where a value of the power density E was calculated using the data from the wind sites discussed in Section 2. To obtain this figure, the values of the rated, cut-in and cut-out wind speeds were changed independently to study the individual effect of each parameter on the energy production values.The values presented in this section will be used to define the power curve expression used in the following sections.
In all the methods used in this paper, the lower wind speeds (calms) were treated separately [2,6] to improve the estimation results.Calm wind speeds were removed from the measurements; the fitting method was then applied to the remaining data.Finally, the calms speeds were re-introduced into the results to properly calculate the energy yield and other parameters.

Goodness of Fit
To determine if a theoretical PDF is suitable to describe the wind speed data, several indicators can be used.For each wind site and fitting method, the following indicators were considered: • The relative mean wind speed error (error vm ) and the mean wind speed data (v m ) are compared with the resulting mean wind speed ( m v′ ) from the fitted PDF using: • The relative error of the available power density (error Ew ): where E w is calculated from the measured data; and w E′ is calculated from the estimated PDF.
• The coefficient of determination of the wind speed distribution (R 2 ): where j fr′ is the estimated relative frequency of bin "j"; and 1 1 N j j fr fr N = =  is the mean of the f rj values.In this case, the goodness of fit is better when the coefficient approaches 1.When R 2 is applied to the distribution of the wind power density, the following indicator is obtained: where w e is the mean of the e w (v j ) values.
• Root mean square error (RMSE): ( ) The goodness of fit is better when the RMSE approaches zero.Goodness of fit parameters related to hypothesis testing methods, such as the Chi-square or Anderson-Darling methods, were not used in this work because their values strongly depend on the number of data points, which makes it difficult to compare results from different wind sites.

Proposed Indicators of Fitness
The energy produced by WTGs should be taken into account when Weibull PDFs are used in wind energy applications.For this reason, a set of indicators are proposed in the following paragraphs to include the energy output of WTGs in indicator calculations, which is done to increase the independence of the results from those obtained supposing a particular WTG or a reduced set [20,25,26].
For the aforementioned purpose, a set of power curves is considered, which is defined by selecting different rated wind speeds using Equation (14).In this case, the rated wind speeds v r is defined between 10 m/s and 17 m/s, according to the values shown in Section 5.4 (see Figure 4) and can be represented as: where N r is the number of rated wind speeds used in the histogram.Therefore, at each wind site, a set of wind energy distributions is calculated using the different power curves defined by the rated wind speeds in Equation (22).The wind energy distributions obtained from the data and from the fitted PDF are compared using the proposed indicators: • Relative error of power density (error mE ), which is obtained as the mean value of power density errors at different rated wind speeds using the wind speed data and fitted PDF: where l E and l E′ are the total power density values calculated at a specific wind site using the wind speed data and the estimated PDF with the rated wind speed v rl , respectively.
• Coefficient of determination of the power density distribution ( 2 me R ), which is obtained as the mean value of the R 2 values calculated with the power density distribution at different rated wind speeds using the wind speed data and fitted PDF: where e l and l e′ are the power density distributions calculated at a specific wind site using the data and the estimated PDF with the rated wind speed v rl , respectively, and l e is the mean for the e l (v j ) values.

Proposed Estimation Method: Part Density Energy Method (PDEM)
Regarding the methods to estimate the Weibull parameters, despite the fact that the estimation of energy produced by WTGs is one of the primary objectives in evaluating wind sites, the only method that partially takes it into account is the power density method (PDM), which uses the available energy associated with a wind speed distribution.In this context, the method, named the PDEM, which considers the typical behavior of power curves, is proposed in this paper, which is accomplished using a method that minimizes the following index: where s is the class represented by the wind speed v s , which is selected as the limit to calculate the energy at low and high wind speeds.Therefore, using this function, the energy produced at low wind speeds, where WTGs typically follow a cubic equation [18,32], and the energy produced at high wind speeds, where WTGs typically follow a flat power curve, are considered separately.The index J is minimized using the Nelder-Mead simplex method [40].
The selection of v s in Equation ( 25) was done using the results shown in Figure 6.In this figure, the mean value and the STD of the values of error mE obtained from the PDEM method at all wind sites are represented against different v s values.As can be seen, at the selected v s , the mean error is at its minimum value.Therefore, a value of v s equal to 12 m/s was selected.

Parameter Estimation Results
After choosing all the estimation methods (MLE, MMLE, LSQM, MM, PDM and PDEM), the methods were used to estimate the Weibull PDF parameters for each wind site.An example of the estimated PDF and distribution of the wind power density is plotted in Figures 7 and 8.

PDEM
The Weibull parameters obtained for each wind site using the different fitting methods are shown in Table 2 and in Figure 9. Once the parameters of the Weibull distribution are obtained, the goodness of fit indicators is calculated.The results for each wind site are displayed in Figures 10-16.A summary of the results is presented in Table 3, which lists the mean and standard deviation (STD) of the indicator of the fitness values obtained for all wind sites for each estimation method.
To evaluate the fitting methods, all the fitness indicators should be taken into account because, for example, higher R 2 values, when fitting wind, wind energy or energy, do not imply lower wind, wind energy or energy errors [26].The following are the main conclusions regarding the fitting methods: • The proposed PDEM method exhibited the best behavior when the R 2 values for wind power density and power density distributions were considered.In addition, this method has an acceptable behavior when the relative error of the wind power density and power density were taken into account.The overall behavior of the proposed method is extremely satisfactory.
• MM exhibited the best behavior in wind distribution fitting according to the R 2 and RMSE values [19].Furthermore, the method resulted in acceptable values for all the indicators considered in this paper.
• PDM's behavior is satisfactory with respect to estimating the mean wind speed; however, the method failed when the wind power production was considered [19,24].
• LSQM, in terms of relative error, strongly depended on the wind site data, as can be seen in the high values shown by the STD of its errors (error vm , error Ew and error mE ) [19,20,22].
In conclusion, the proposed PDEM is the most suitable method when the focus of the analysis is on the energy produced by WTGs.Nevertheless, other methods, particularly MM, exhibited satisfactory results in terms of energy fitness.
To evaluate the robustness of the proposed PDEM, the behaviour of the different fitting methods with an estimated wind speed at different hub heights is shown in Appendix B, where it can be seen that similar conclusions are achieved.

Conclusions
This paper presents an analysis of wind speed data based on using fitting curve methods to obtain the parameters of the Weibull PDF.The most widely used methods were selected for this analysis.Furthermore, a method, called PDEM, which takes into account the power density distribution and the typical performance of WTGs, is presented.
The results of the fitting methods in obtaining the Weibull parameters were evaluated by a set of indicators defined from wind speeds and wind power density distributions.Additionally, the indicators, error mE and 2 me R , which consider the behavior of WTGs in terms of energy production, are introduced.
As the primary conclusion, the proposed PDEM method exhibited the best results when the energy produced by the WTGs is considered.Furthermore, its result when all indicators are taken into account, are extremely satisfactory.Nevertheless, other methods, particularly MM, exhibited satisfactory results in terms of energy fitness.In this paper, wind speed data from weather stations from a specific region, northwest Spain (Galicia), were used.
( ) ( ) where v i are the wind speed data values.Using data from a histogram, the modified maximum likelihood method (MMLM) results [20] in the following: ( ) ( ) where v j are the wind speed representatives of each bin of the histogram.

A2. LSQM
The LSQM, also known as the Weibull plot, is based on logarithmic transformations applied to the Weibull cumulative distribution function F(v) and thus, can be represented by a straight line: where the straight line can be written as: ( ) The Weibull parameters can be obtained from: Finally, the line parameters are obtained from: where x m and y m are the mean values of x and y, respectively; and x i and y i are the values obtained using Equation (A5) with each wind speed data set.The LSQM is extremely popular due to its simplicity; however, the logarithmic transformations used during the calculation tend to cause some inaccuracy [20].

A3. MM
The MM is based on obtaining the Weibull parameters using certain statistical moments calculated using wind speed data.When the mean wind speed (v m ) and the STD (σ) of wind data are used, the following relationship can be obtained [24]: The shape parameter k can be calculated from this equation by the Newton-Raphson (NR) method: where z = 1/k.Nevertheless, the following approximation typically gives acceptable results: ( ) The comparison between the approximate solution using Equation (A10) and the solution solving Equation (A9) by the Newton Raphson method is shown in Figure 17.To avoid errors, the approximated solution was not used in this paper.

A4. PDM
The PDM uses the Energy Pattern Factor [24,25]: This factor relates to the shape parameter k by means of: ( ) ( ) Using the NR method to solve this equation, the following expressions are obtained: The following approximated expression can be used, which assumes that E pf is typically between 1.45 and 4.4: Both results from using Equations (A13) and (A14) are shown in Figure 17.To avoid errors, the approximated solution was not used in this paper.

B. Behavior of the Proposed PDEM Method at Different Hub Heights
In order to evaluate the robustness of the proposed PDEM method, it has been evaluated using wind speeds at different hub heights.For this purpose, wind speeds have been estimated using a logaritmic wind profile with a roughness length of 0.05 m, which is a common value in wind farms [9,12].The resulting mean values of the indicators of the fitness (error vm , error Ew , error mE , R 2 ew , R 2 me and R 2 ) for heights between 10 m and 150 m can be seen in Figure 18.As a conclusion, the overall relative behaviour of the proposed PDEM method does not significantly change with the hub height.

Figure 1 .
Figure 1.Meteogalicia weather stations, where the size of each circle is proportional to the measured mean wind speed, and the number inside each circle is the site number.The circles at the top left portion of the figure indicate the scale.

Figure 3 .
Figure 3. Power density E (at rated speeds v r of 11 m/s and 15 m/s) and wind power density E w for different Weibull parameters.

Figure 4 .
Figure 4. Distribution of the cut-in, rated and cut-out wind speeds for wind turbine generators (WTGs).
wind speed in m/s rel.freq.%

Figure 5 .
Figure 5. Candlestick chart with the power density at the considered wind sites for different values of rated, cut-in and cut-out wind speeds.

Figure 6 .
Figure 6.Evolution of the mean and standard deviation of the relative error error mE at different selected wind speeds v s .

Figure 7 .
Figure 7. Histogram derived from the estimated Weibull probability density function (PDF) compared with the histogram of the wind speed data (bar diagram) at wind site n° 21(O Cebreiro).MLE: maximum likelihood method; MMLE: modified maximum likelihood estimation; LSQM: least square method; MM: moment method; PDM: power density method; and PDEM: part density energy method.

Figure 8 .
Figure 8. Distribution of the wind power density derived from the estimated Weibull PDF compared with that derived from the wind speed data (bar diagram) at wind site n° 21 (O Cebreiro).

Figure 9 .
Figure 9. Representation of the mean wind speed calculated from data, scale and shape parameters for all wind sites using the different estimation methods.

Figure 10 .
Figure 10.Representation of the error vm at the different wind sites.

Figure 11 .
Figure 11.Representation of the error Ew at the different wind sites.

Figure 12 .
Figure 12.Representation of R 2 at the different wind sites.

Figure 13 .Figure 14 .
Figure 13.Representation of 2 ew R at the different wind sites

Figure 15 .Figure 16 .
Figure 15.Representation of each component of error mE at different rated wind speeds at the different wind sites.For each wind site, the final value of error mE is the mean value of all the components.

Figure 17 .
Figure 17.Shape parameter estimated with MM and density power method (DPM).

Figure 18 .
Figure 18.Mean values of the indicators of the fitness (error vm , error Ew , error mE , R 2 ew , R 2 me and R 2 ) at different hub heights.

Table 1 .
Meteorological stations (No.: wind site number; MWS: mean wind speed in m/s; and ND: amount of data in years).

Table 2 .
Weibull PDF parameters: scale parameter c in m/s and unitless shape parameter k.

Table 3 .
Summary of the results (best values of the mean and standard deviation (STD) for each fitness indicator are displayed in bold).