jmquintana79
9/5/2017 - 6:37 AM

## Wind Analysis Plot

Plot wind rose and wind distribution. Furthermore, return a Wind Speed Frequency Matrix.

``````"""
WIND ANALYSIS:

Plot the wind rose and wind distribution. It is required a DF with hwd and hws data.
"""
def analysis_wind(DATA,svar_hws,svar_hwd,stitle=''):

## CALCULATE SECTORS
def setsector(vdir):
for k, v in dsector.items():
if k==0:
if (vdir>=v[1]) or (vdir<v[2]): return v[0]
else:
if (vdir>=v[1]) and (vdir<v[2]): return v[0]

## SECTOR DEFINITION (16 sectors)

# initialize
dsector = dict(); disector = dict()
lsectors = ['N','NNE','NE','ENE','E','ESE','SE','SSE','S','SSW','SW','WSW','W','WNW','NW','NNW']
llimits =[348.75,11.25],[11.25,33.75],[33.75,56.25],[56.25,78.75],[78.75,101.25],[101.25,123.75],[123.75,146.25],[146.25,168.75],[168.75,191.25],[191.25,213.75],[213.75,236.25],[236.25,258.75],[258.75,281.25],[281.25,303.75],[303.75,326.25],[326.25,348.75]
# loop of angles
for i,isector in enumerate(lsectors[:]):
# build dictionary
dsector[i] = [isector,llimits[i][0],llimits[i][1]]
disector[isector] = i

# collect sector
lsector = [v[0] for k,v in dsector.items()]
# calculate sector
DATA['sector'] = DATA[svar_hwd].apply(lambda col: setsector(col))

## BIN DEFINITION

# set bin
def setbin(vspeed):
for k, v in dbin.items():
if (vspeed>=v[1]) and (vspeed<v[2]): return int(k)

ibin = 0
dbin = dict()
for il in range(16):
# set limits
il1 = il
if il==15: il2 = 9999
else: il2 = il+1
# define bin and store
dbin[ibin]=['%s<=hws<%s'%(il1,il2),il1,il2]
# define number of bin
ibin += 1

# calculate bin
DATA['bin'] = DATA[svar_hws].apply(lambda col: setbin(col))

""" CREATE HWS FREQUENCY MATRIX """

# initialize
lbin = list(range(16))
ntotal = len(DATA)
dfreq = dict(); davg = dict()
# loop of sectors
for isector in lsector:
# initialize
lfreq = list(); lavg = list()
# loop of bins
for ibin in lbin:
# collect data
AUX = DATA[(DATA['sector']==isector) & (DATA['bin']==ibin)]
# calculate frequency (%)
ifreq = len(AUX) * 100. / ntotal
# calculate avg of hws
iavg = np.mean(AUX[svar_hws].values)
# store
lfreq.append(ifreq); lavg.append(iavg)
# clean
del(AUX)
# store by sector
dfreq[isector] = lfreq
davg[isector] = [np.nanmean(lavg)]

## FINAL DF CREATION

# create frequency matrix
WSFREQ = pd.DataFrame.from_dict(dfreq, orient='columns')
WSFREQ['index'] = [v[0] for k,v in dbin.items()]
WSFREQ = WSFREQ.set_index('index')
WSFREQ['ALL'] = WSFREQ.sum(axis=1)
# create avg hws
WSAVG = pd.DataFrame.from_dict(davg, orient='columns')
WSAVG['ALL'] = WSAVG.mean(axis=1)
WSAVG['index'] = ['avg_hws(m/s)']
WSAVG = WSAVG.set_index('index')
# create sum freq
WSFREQSUM = pd.DataFrame(WSFREQ.sum(axis=0)).transpose()
WSFREQSUM['index'] = ['sum_freq(%)']
WSFREQSUM = WSFREQSUM.set_index('index')
# append total
WSFREQ = WSFREQ.append(WSAVG)
WSFREQ = WSFREQ.append(WSFREQSUM)
# clean
del(WSAVG); del(WSFREQSUM)

""" PLOT WIND ROSE """

import matplotlib.pyplot as plt
from matplotlib import rcParams
from matplotlib.colors import ListedColormap
import matplotlib as mpl

# CREATE A TWIN POLAR AXIS
def polar_twin(ax):
label='twin', frameon=False,
theta_direction=ax.get_theta_direction(),
theta_offset=ax.get_theta_offset())
ax2.xaxis.set_visible(False)

# There should be a method for this, but there isn't... Pull request?
ax2._r_label_position._t = (22.5 + 180., 0.0)
ax2._r_label_position.invalidate()

# Bit of a hack to ensure that the original axes tick labels are on top of
# whatever is plotted in the twinned axes. Tick labels will be drawn twice.
for label in ax.get_yticklabels():
ax.figure.texts.append(label)

return ax2

## COLLECT DATA
PLOT = WSFREQ.transpose()[['avg_hws(m/s)','sum_freq(%)']][:-1]
PLOT.columns = ['avg','freq']
PLOT['isector'] = [disector[i] for i in PLOT.index.values]
PLOT.sort_values(['isector'], ascending=[1], inplace=True)

## DATA PREPARATION

# collect total data
values = list(PLOT.avg.values)
#freq = list(PLOT.freq.values)
lsector = list(PLOT.index.values)
# clean
del(PLOT)

# calculate total freq per sector for 4 defined bins
# bin1: [0,4)  bin2: [4,8)  bin3: [8,12)  bin4: [12,inf)
ib1 = (0,3); ib2 = (4,7); ib3 = (8,11); ib4 = (12,15)
lfsector1 = list(); lfsector2 = list(); lfsector3 = list(); lfsector4 = list()
for isector in lsector:
lfsector1.append(np.nansum(WSFREQ[ib1[0]:ib1[1]+1][isector]))
lfsector2.append(np.nansum(WSFREQ[ib2[0]:ib2[1]+1][isector]))
lfsector3.append(np.nansum(WSFREQ[ib3[0]:ib3[1]+1][isector]))
lfsector4.append(np.nansum(WSFREQ[ib4[0]:ib4[1]+1][isector]))

# calculate ylim
freq = np.array(lfsector1)+np.array(lfsector2)+np.array(lfsector3)+np.array(lfsector4)
ylim = int(np.max(freq))+2

# define angles of sectors for labels of rose chart
langle = list(np.arange(0,360,22.5))

# define x
numpoints = len(lsector)
theta = [n / float(numpoints) * 2 * pi for n in range(numpoints)]

# Because our chart will be circular we need to append a copy of the first
# value of each list at the end of each list with data
values += values[:1]
theta = list(theta)
thetaf = copy.deepcopy(theta)
theta += theta[:1]

## create figure object
params = dict(projection='polar', theta_direction=-1, theta_offset=np.pi/2)
fig, ax = plt.subplots(subplot_kw=params,figsize=(6,7))

## PLOT CHART OF FREQUENCIES

# plot
width = 0.25
ax.bar(thetaf,lfsector1,width, color='blue', alpha=0.8)
ax.bar(thetaf,lfsector2,width,bottom=np.array(lfsector1), color='cyan', alpha=0.8)
ax.bar(thetaf,lfsector3,width,bottom=np.array(lfsector1)+np.array(lfsector2), color='orange', alpha=0.8)
ax.bar(thetaf,lfsector4,width,bottom=np.array(lfsector1)+np.array(lfsector2)+np.array(lfsector3), color='red', alpha=0.8)
# define limits
ax.set_ylim([0, ylim])
# set axi tick labels
plt.setp(ax.get_yticklabels(), color='blue')
plt.yticks(range(0,ylim+4,4),fontsize=10,rotation='horizontal')
ax.set_yticklabels(list(range(0,ylim+4,4)))

## PLOT CHART OF Avg. WS

# define second axis as twin
ax2 = polar_twin(ax)
# plot
ax2.plot(theta, values, color='black', alpha=1.)
# set limits
ax2.set_ylim([0, ylim])
# set axi tick labels
plt.setp(ax2.get_yticklabels(), color='red')
plt.yticks(range(0,ylim+4,4),fontsize=10,rotation='horizontal')
ax2.set_yticklabels(list(range(0,ylim+4,4)))

# set axis labels
ax2.text(radians(22.5+7.5), ylim + 5, '(m/s)', size=10,color='blue', horizontalalignment="center", verticalalignment="center")
ax.text(radians(202.5-7.5), ylim -0, '(%)', size=10,color='red', horizontalalignment="center", verticalalignment="center")

# set ytick labels of rose chart
ax.set_thetagrids(langle,labels=lsector,fontsize=16, weight="normal", color="black")

## LEGEND

# add new axis for legend object
ax3 = fig.add_axes([0.5-0.25/2., 0.08, 0.25, 0.025])

# value colors.
cmap = mpl.colors.ListedColormap(['blue', 'cyan', 'orange','red'])
cmap.set_over('0.25')
# create boundary object
bounds = [0, 4, 8, 12, 15]
norm = mpl.colors.BoundaryNorm(bounds, cmap.N)
# create color bar
cb = mpl.colorbar.ColorbarBase(ax3, cmap=cmap,
norm=norm,
boundaries=bounds + [13],
extend='max',
ticks=bounds,
spacing='proportional',
orientation='horizontal')
# set label of color bar

# title
if len(stitle)>0: title = ': %s'%stitle
else: title = stitle
ax.set_title('Wind Rose%s'%title,fontsize=15,y=1.1)

# display
plt.show()

""" PLOT BAR FREQ """

import matplotlib.pyplot as plt

## COLLECT DATA
BAR = WSFREQ[['ALL']].iloc[0:16]

## PREPARE DATA
nrow = len(BAR)
x = np.array(list(range(nrow)))+0.5
y = BAR.ALL.values
ylim = int(np.max(y))+2
# clean
del(BAR)

## PLOT

# create objects
fig, ax = plt.subplots(figsize=(6,7))

# plot bars
width = 1
ax.bar(x,y,width, color='r', edgecolor='black')

# limits of axes
ax.set_xlim([0,16])
ax.set_ylim([0,ylim])

# resolution of axis ticks
plt.xticks(list(range(nrow)),fontsize='small',rotation='horizontal')
ax.set_xticklabels(list(range(nrow)))
plt.yticks(range(0,ylim+4,4),fontsize='small',rotation='horizontal')
ax.set_yticklabels(list(range(0,ylim+4,4)))

# set axis labels
plt.xlabel("WS(m/s)",fontsize=10)
ax.set_ylabel("freq(%)",fontsize=10)

# title
if len(stitle)>0: title = ': %s'%stitle
else: title = stitle
plt.title('Wind Speed Distribution%s'%title,fontsize=16)

# set labels on points for scatter chart
lvalue = ['%.1f'%iy for iy in y]
for i, svalue in enumerate(lvalue):
if float(svalue)>=10.: ax.annotate(svalue, (x[i]-0.48,y[i]+0.25),fontsize=9,color='red',rotation='horizontal')
else: ax.annotate(svalue, (x[i]-0.3,y[i]+0.25),fontsize=9,color='red',rotation='horizontal')

# display
plt.show()

# close plot
plt.cla()   # Clear axis
plt.clf()   # Clear figure
plt.close() # Close a figure window

# return data matrix
return WSFREQ``````