Analysing wind data with R

Question

Hi i am analaysing wind data for estimating energy from a wind turbine.
I have taken 10 years of wind data and graphed a histogram;
my second stage was to fit a Weibull distribution to the data.
I used R with the package lmom to compute the Weibul shape and scale this is the code i used:

>library(lmom)    
wind.moments<-samlmu(as.numeric(pp$WS))      
moments<-pelwei(wind.moments)     
x.wei<-rweibull(n=length(pp$WS), shape=moments["delta"], scale=moments["beta"])    
hist(as.numeric(pp$WS), freq=FALSE)    
lines(density(x.wei), col="red", lwd=4)    

It seems like there is some lag between the data and the density function; can you help me with this? Another question is can you help me in calculating the anual energy from the density function?

thank you

Answer 1

lmom function pelwei fits a three parameter Weibull distribution, with location, scale and shape parameters. rweibull generates random numbers for a two-parameter Weibull distribution. You need to subtract the location parameter moments["zeta"]. That should give a better fit, but it doesn't appear it will give a good fit to your particular data.

I notice http://www.reuk.co.uk/Wind-Speed-Distribution-Weibull.htm says "Wind speeds in most of the world can be modelled using the Weibull Distribution.". Perhaps you're just unlucky and live in a part of the world where they can't!

Answer 2

I recreated your plot with data from http://hawaii.gov/dbedt/ert/winddata/krab0192.txt (I took the 1200 measurements). I got a decent fit of the data, generally using your code:

library(lmom)

daten <- read.delim("wind.txt")
wind.avg <- na.omit(as.numeric(daten[,"X12"]))
wind.moments<-samlmu(wind.avg)
moments<-pelwei(wind.moments)
x.wei<-rweibull(n=length(wind.avg), shape=moments["delta"], scale=moments["beta"])
hist(as.numeric(wind.avg), freq=FALSE)
lines(density(x.wei), col="red", lwd=4)

Sorry, I'm not shure were your problem could be, but I think you should be able to fit weibull to your data. What makes me suspicious is the bell-curve of your density plot, I have no idea where that came from.

Here are the moments I generated:

wind.moments

       l_1         l_2         t_3         t_4 
15.17287544  4.80372580  0.14963501  0.06954438

moments

     zeta      beta     delta 
 0.516201 16.454233  1.745413 

WTR to the annual output: I suppose I'd generate discrete values for the probability density function, multiply these values with the output function and sum it up. Alternatively, you could just use your raw data, multiply the values with the output function, sum it up and calculate the annual average, you should control for seasonality in a suitable way (e. g. make sure to use whole years, or to weight accordingly).

Here is the uncontrolled output (using the formula from http://www.articlesbase.com/diy-articles/determining-wind-turbine-annual-power-output-a-simple-formula-based-upon-blade-diameter-and-average-wind-speed-at-your-location-513080.html)

years  <- length(wind.avg)/365
diameter <- 150
Power = (0.01328*diameter^2)*((wind.avg)^3)
(annual.power <- sum(Power)/years)
[1] 791828306

Answer 3

Here's a recent post at SO on wind turbines. My answer on that link has three links that you might be interested in:

http://stackoverflow.com/questions/4843194/r-language-sorting-data-into-ranges-averaging-ignore-outliers/4844783#4844783

I just checked one of the Weibull links in the above SO answer. For some reason, the link is down. Here are some links that provide the same basic information:

http://www.gso.uri.edu/ozone/

http://www.weru.ksu.edu/new_weru/publications/pdf/Comparison%20of%20the%20Weibull%20model%20with%20measured%20wind%20speed%20distributions%20for%20stochastic%20wind%20genera.pdf

http://www.kfupm.edu.sa/ri/publication/cer/41_JP_Weibull_parameters_for_wind_speed_distribution_in_Saudi_Arabia.pdf

http://journal.dogus.edu.tr/13026739/2008/cilt9/sayi1/M00195.pdf

http://www.eurojournals.com/ejsr_26_1_01.pdf

Also, from the power generated from wind, the seasonality is obvious.

Answer 4

I'm not sure if somebody has already made this point, but pelwei can actually be forced to work as a 2 parameter weibull function by adding in a fixed bound.

Insead of calling moments<-pelwei(wind.moments) you should simply call moments<-pelwei(wind.moments,bound=0)

you can always check what the zeta value is. If it's not 0 and you're using dweibull, you need to do something about it.