Code covered by the BSD License  

Highlights from
Example files for "Programming with MATLAB" Webinar

image thumbnail

Example files for "Programming with MATLAB" Webinar


Jiro Doke (view profile)


Example files for "Programming with MATLAB" Webinar, first delivered October 3, 2013.

Wind Turbine Data Analysis

Wind Turbine Data Analysis

This demo analyzes wind data measured on a meteorological observation tower to see if the location is a good prospect for a wind turbine. Data is from three different wind sensors at 80m. Temperature is also recorded at 3m height. Data is logged every hour.


Read in Turbine Data from Text File

% Function autogenerated from import tool
filename = 'winddata.txt';

Visualize Wind Speed and Temperature

figure; plot(time, [velS1, velS2, velS3])
datetick('x', 'mmm-yy')
ylabel('wind velocity at 80m (m/s)')
title(['Wind Data from ', filename])

figure; plot(time, tempS, 'r')
datetick('x', 'mmm-yy')
ylabel('temperature (\circC)')
title(['Temperature Data from ', filename])

Average Wind Speed for Different Sensors

velAvg = mean([velS1, velS2, velS3], 2);

Determine Icing Conditions

Remove any readings effected by icing which results in extremely inaccurate values biased towards zero. Icing conditions are when tavg < tIce and vavg < vIce.

tIce = 2;
vIce = 1;

% Comparing to the critical values
idxTIce = tempS < tIce;
idxVIce = velAvg < vIce;

idxIce = idxTIce & idxVIce ;

% Plotting when icing occurs
figure; plot(time, tempS);
hold on
plot(time(idxIce), tempS(idxIce), 'r*');
datetick('x', 'mmm-yy')
ylabel('temperature (\circC)')
title (['Removed Readings from ', filename, ' due to Icing'])

% Remove values related to icing
time(idxIce)   = [];
tempS(idxIce)  = [];
velAvg(idxIce) = [];

Distribution of Wind Speeds at Hub Height

Use a weibull distribution to fit the distribution of wind speeds which is known to often give a good fit to wind speed distributions.

% Plotting a histogram of wind speeds
dv = 0.5;
vbins = 0:dv:ceil(max(velAvg));

% Plotting a histogram of wind speeds
figure; hist(velAvg, vbins);
xlabel('wind velocity (m/s)'); ylabel('count');

% Calculate probability of wind speed range
nelements = hist(velAvg, vbins);
probvbins = nelements/sum(nelements);        % Probability of a given velocity range

Defining the Turbine Power Curve

To calculate average turbine power and capacity factor, we need to make some assumptions regarding the wind turbine model and its power curve. We will assume a 1MW wind turbine and the following turbine power curve.

% Turbine power curve coefficents
prated = 1e6;     % wind turbine rated power (W)
vin = 2;          % cut-in speed (m/s)
vr = 14;          % rated output speed (m/s)
vout = 25;        % cut-out speed (m/s)

% Calculating power curve
powervbins = prated*(vbins.^2 - vin^2)/(vr^2 - vin^2);
powervbins(vbins <= vin) = 0;
powervbins(vbins > vout) = 0;
powervbins(vbins >= vr & vbins <= vout) = prated;

% Plotting power curve
figure; plot(vbins, powervbins, '*')
xlabel('velocity (m/s)');
ylabel('turbine power (W)');
title('Turbine Power Curve')

Calculating Average Turbine Power and Capacity Factor

Capacity factor is ratio of the actual output of a turbine over a period of time and its potential output if it had operated at full capacity the entire time. Typical capacity factor range from 20-–50% depending on location and wind turbine.

% Integrate power at given velocity * velocity probability distribution
% function over range of possible velocities

% Calculate average power by summing the products of vel probability and power
avgPower = sum(probvbins .* powervbins);    % (W)

% Calculating capacity factor (average power / rated power)
cf = avgPower / prated;

disp(['Assumed wind turbine rated power (MW): ', num2str(prated/1e6)]);
disp(['Averaged turbine power (kW): ', num2str(avgPower/1e3)]);
disp(['Capacity factor (%): ', num2str(cf*100)]);
Assumed wind turbine rated power (MW): 1
Averaged turbine power (kW): 270.5152
Capacity factor (%): 27.0515

Contact us