diff --git a/wetb/fatigue_tools/fatigue.py b/wetb/fatigue_tools/fatigue.py
index e26089ce4f5c6717651908ac1e35ef918d52fb9b..bde2d087b59a418c8fa6040f5a3c126e9330277b 100644
--- a/wetb/fatigue_tools/fatigue.py
+++ b/wetb/fatigue_tools/fatigue.py
@@ -70,7 +70,6 @@ def eq_load(signals, no_bins=46, m=[3, 4, 6, 8, 10, 12], neq=1, rainflow_func=ra
return [[np.nan] * len(np.atleast_1d(m))] * len(np.atleast_1d(neq))
-
def eq_load_and_cycles(signals, no_bins=46, m=[3, 4, 6, 8, 10, 12], neq=[10 ** 6, 10 ** 7, 10 ** 8], rainflow_func=rainflow_windap):
"""Calculate combined fatigue equivalent load
@@ -109,7 +108,6 @@ def eq_load_and_cycles(signals, no_bins=46, m=[3, 4, 6, 8, 10, 12], neq=[10 ** 6
return eq_loads, cycles, ampl_bin_mean, ampl_bin_edges
-
def cycle_matrix(signals, ampl_bins=10, mean_bins=10, rainflow_func=rainflow_windap):
"""Markow load cycle matrix
@@ -132,7 +130,8 @@ def cycle_matrix(signals, ampl_bins=10, mean_bins=10, rainflow_func=rainflow_win
Returns
-------
cycles : ndarray, shape(ampl_bins, mean_bins)
- A bi-dimensional histogram of load cycles(full cycles). Amplitudes are histogrammed along the first dimension and mean values are histogrammed along the second dimension.
+ A bi-dimensional histogram of load cycles(full cycles). Amplitudes are\
+ histogrammed along the first dimension and mean values are histogrammed along the second dimension.
ampl_bin_mean : ndarray, shape(ampl_bins,)
The average cycle amplitude of the bins
ampl_edges : ndarray, shape(ampl_bins+1,)
diff --git a/wetb/fatigue_tools/rainflowcounting/rfc_hist.py b/wetb/fatigue_tools/rainflowcounting/rfc_hist.py
deleted file mode 100644
index 65d3c00945c8f97f731af59972c4a8b2f55a7b25..0000000000000000000000000000000000000000
--- a/wetb/fatigue_tools/rainflowcounting/rfc_hist.py
+++ /dev/null
@@ -1,53 +0,0 @@
-from __future__ import division
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import absolute_import
-from future import standard_library
-standard_library.install_aliases()
-import numpy as np
-def rfc_hist(sig_rf, nrbins=46):
- """Histogram of rainflow counted cycles
-
- hist, bin_edges, bin_avg = rfc_hist(sig, nrbins=46)
-
- Divide the rainflow counted cycles of a signal into equally spaced bins.
-
- Created on Wed Feb 16 16:53:18 2011
- @author: David Verelst
- Modified 10.10.2011 by Mads M Pedersen to elimintate __copy__ and __eq__
-
- Parameters
- ----------
- sig_rf : array-like
- As output by rfc_astm or rainflow
-
- nrbins : int, optional
- Divide the rainflow counted amplitudes in a number of equally spaced
- bins.
-
- Returns
- -------
- hist : array-like
- Counted rainflow cycles per bin, has nrbins elements
-
- bin_edges : array-like
- Edges of the bins, has nrbins+1 elements.
-
- bin_avg : array-like
- Average rainflow cycle amplitude per bin, has nrbins elements.
- """
-
- rf_half = sig_rf
-
- # the Matlab approach is to divide into 46 bins
- bin_edges = np.linspace(0, 1, num=nrbins + 1) * rf_half.max()
- hist = np.histogram(rf_half, bins=bin_edges)[0]
- # calculate the average per bin
- hist_sum = np.histogram(rf_half, weights=rf_half, bins=bin_edges)[0]
- # replace zeros with one, to avoid 0/0
- hist_ = hist.copy()
- hist_[(hist == 0).nonzero()] = 1.0
- # since the sum is also 0, the avg remains zero for those whos hist is zero
- bin_avg = hist_sum / hist_
-
- return hist, bin_edges, bin_avg
diff --git a/wetb/hawc2/Hawc2io.py b/wetb/hawc2/Hawc2io.py
index 2ec8deff39bf036f1879869ac3121dce1b5ecbb6..ac8d28c8d43327f10373873b4521febfaf041434 100644
--- a/wetb/hawc2/Hawc2io.py
+++ b/wetb/hawc2/Hawc2io.py
@@ -90,8 +90,8 @@ class ReadHawc2(object):
Name = []; Unit = []; Description = [];
for i in range(0, self.NrCh):
temp = str(Lines[i + 12][12:43]); Name.append(temp.strip())
- temp = str(Lines[i + 12][43:48]); Unit.append(temp.strip())
- temp = str(Lines[i + 12][49:]); Description.append(temp.strip())
+ temp = str(Lines[i + 12][43:54]); Unit.append(temp.strip())
+ temp = str(Lines[i + 12][54:-1]); Description.append(temp.strip())
self.ChInfo = [Name, Unit, Description]
# if binary file format, scaling factors are read
if Format.lower() == 'binary':
diff --git a/wetb/prepost/windIO.py b/wetb/prepost/windIO.py
index 5c38c195d772b9eb74eba2f86c53c504060572ae..d36d3f3edcdca1983128dd3e1a2ec02d794069be 100755
--- a/wetb/prepost/windIO.py
+++ b/wetb/prepost/windIO.py
@@ -17,10 +17,6 @@ from future import standard_library
standard_library.install_aliases()
from builtins import object
- # always devide as floats
-
-#print(*objects, sep=' ', end='\n', file=sys.stdout)
-
__author__ = 'David Verelst'
__license__ = 'GPL'
__version__ = '0.5'
@@ -38,18 +34,16 @@ import array
import numpy as np
import pandas as pd
-#import sympy
-
# misc is part of prepost, which is available on the dtu wind gitlab server:
# https://gitlab.windenergy.dtu.dk/dave/prepost
from wetb.prepost import misc
# wind energy python toolbox, available on the dtu wind redmine server:
# http://vind-redmine.win.dtu.dk/projects/pythontoolbox/repository/show/fatigue_tools
-from wetb.fatigue_tools.rainflowcounting.rainflowcount import rainflow_astm as rainflow_astm
-from wetb.fatigue_tools.rainflowcounting.rfc_hist import rfc_hist as rfc_hist
+from wetb.hawc2.Hawc2io import ReadHawc2
+from wetb.fatigue_tools.fatigue import eq_load
-class LoadResults(object):
+class LoadResults(ReadHawc2):
"""Read a HAWC2 result data file
Usage:
@@ -126,215 +120,34 @@ class LoadResults(object):
self.file_path = file_path
# remove .log, .dat, .sel extensions who might be accedental left
- if file_name[-4:] in ['.htc','.sel','.dat','.log']:
+ if file_name[-4:] in ['.htc', '.sel', '.dat', '.log']:
file_name = file_name[:-4]
-
+ # FIXME: since HAWC2 will always have lower case output files, convert
+ # any wrongly used upper case letters to lower case here
self.file_name = file_name
- self.read_sel()
- # create for any supported channel the
- # continue if the file has been succesfully read
- if self.error_msg == 'none':
- # load the channel id's and scale factors
- self.scale_factors = self.data_sel()
- # with the sel file loaded, we have all the channel names to
- # squeeze into a more consistant naming scheme
- self._unified_channel_names()
- # only read when asked for
- if readdata:
- # if there is sel file but it is empty or whatever else
- # FilType will not exists
- try:
- # read the binary file
- if self.FileType == 'BINARY':
- self.read_bin(self.scale_factors, usecols=usecols)
- # read the ASCII file
- elif self.FileType == 'ASCII':
- self.read_ascii(usecols=usecols)
- else:
- print('='*79)
- print('unknown file type: ' + self.FileType)
- print('='*79)
- self.error_msg = 'error: unknown file type'
- self.sig = []
- except:
- print('='*79)
- print('couldn\'t determine FileType')
- print('='*79)
- self.error_msg = 'error: no file type'
- self.sig = []
+ FileName = os.path.join(self.file_path, self.file_name)
- if self.debug:
- stop = time() - start
- print('time to load HAWC2 file:', stop, 's')
+ ReadOnly = 0 if readdata else 1
+ super(LoadResults, self).__init__(FileName, ReadOnly=ReadOnly)
+ self.FileType = self.FileFormat[6:]
+ self.N = int(self.NrSc)
+ self.Nch = int(self.NrCh)
+ self.ch_details = np.ndarray(shape=(self.Nch, 3), dtype='<U100')
+ for ic in range(self.Nch):
+ self.ch_details[ic, 0] = self.ChInfo[0][ic]
+ self.ch_details[ic, 1] = self.ChInfo[1][ic]
+ self.ch_details[ic, 2] = self.ChInfo[2][ic]
- def read_sel(self):
- # anticipate error on file reading
- try:
- # open file, read and close
- go_sel = os.path.join(self.file_path, self.file_name + '.sel')
- FILE = opent(go_sel, "r")
- self.lines = FILE.readlines()
- FILE.close()
- self.error_msg = 'none'
-
- # error message if the file does not exists
- except:
- # print(26*' ' + 'ERROR'
- print(50*'=')
- print(self.file_path)
- print(self.file_name + '.sel could not be found')
- print(50*'=')
- self.error_msg = 'error: file not found'
-
- def data_sel(self):
-
- # scan through all the lines in the file
- line_nr = 1
- # channel counter for ch_details
- ch = 0
- for line in self.lines:
- # on line 9 we can read following paramaters:
- if line_nr == 9:
- # remove the end of line character
- line = line.replace('\n','').replace('\r', '')
-
- settings = line.split(' ')
- # delete all empty string values
- for k in range(settings.count('')):
- settings.remove('')
-
- # and assign proper values with correct data type
- self.N = int(settings[0])
- self.Nch = int(settings[1])
- self.Time = float(settings[2])
- self.FileType = settings[3]
- self.Freq = self.N/self.Time
-
- # prepare list variables
- self.ch_details = np.ndarray(shape=(self.Nch,3),dtype='<U100')
- # it seems that float64 reeds the data correctly from the file
- scale_factors = scipy.zeros(self.Nch, dtype='Float64')
- #self.scale_factors_dec = scipy.zeros(self.Nch, dtype='f8')
- i = 0
-
- # starting from line 13, we have the channels info
- if line_nr > 12:
- # read the signal details
- if line_nr < 13 + self.Nch:
- # remove leading and trailing whitespaces from line parts
- self.ch_details[ch,0] = str(line[12:43]).strip() # chID
- self.ch_details[ch,1] = str(line[43:54]).strip() # chUnits
- self.ch_details[ch,2] = str(line[54:-1]).strip() # chDescr
- ch += 1
- # read the signal scale parameters for binary format
- elif line_nr > 14 + self.Nch:
- scale_factors[i] = line
- # print(scale_factors[i]
- #self.scale_factors_dec[i] = D.Decimal(line)
- i = i + 1
- # stop going through the lines if at the end of the file
- if line_nr == 2*self.Nch + 14:
- self.scale_factors = scale_factors
-
- if self.debug:
- print('N ', self.N)
- print('Nch ', self.Nch)
- print('Time ', self.Time)
- print('FileType', self.FileType)
- print('Freq ', self.Freq)
- print('scale_factors', scale_factors.shape)
-
- return scale_factors
- break
+ ChVec = [] if usecols is None else usecols
- # counting the line numbers
- line_nr = line_nr + 1
+ self._unified_channel_names()
+ if readdata:
+ self.sig = super(LoadResults, self).__call__(ChVec=ChVec)
- def read(self, usecols=False):
- """
- This whole LoadResults needs to be refactered because it is crap.
- Keep the old ones for backwards compatibility
- """
-
- if self.FileType == 'ASCII':
- self.read_ascii(usecols=usecols)
- elif self.FileType == 'BINARY':
- self.read_bin(self.scale_factors, usecols=usecols)
-
- def read_bin(self, scale_factors, usecols=False):
- if not usecols:
- usecols = list(range(0, self.Nch))
- fid = open(os.path.join(self.file_path, self.file_name) + '.dat', 'rb')
- self.sig = np.zeros( (self.N, len(usecols)) )
- for j, i in enumerate(usecols):
- fid.seek(i*self.N*2,0)
- self.sig[:,j] = np.fromfile(fid, 'int16', self.N)*scale_factors[i]
-
- def read_bin_old(self, scale_factors):
- # if there is an error reading the binary file (for instance if empty)
- try:
- # read the binary file
- go_binary = os.path.join(self.file_path, self.file_name) + '.dat'
- FILE = open(go_binary, mode='rb')
-
- # create array, put all the binary elements as one long chain in it
- binvalues = array.array('h')
- binvalues.fromfile(FILE, self.N * self.Nch)
- FILE.close()
- # convert now to a structured numpy array
- # sig = np.array(binvalues, np.float)
-# sig = np.array(binvalues)
- # this is faster! the saved bin values are only of type int16
- sig = np.array(binvalues, dtype='int16')
-
- if self.debug: print(self.N, self.Nch, sig.shape)
-
-# sig = np.reshape(sig, (self.Nch, self.N))
-# # apperently Nch and N had to be reversed to read it correctly
-# # is this because we are reading a Fortran array with Python C
-# # code? so now transpose again so we have sig(time, channel)
-# sig = np.transpose(sig)
-
- # reshape the array to 2D and transpose (Fortran to C array)
- sig = sig.reshape((self.Nch, self.N)).T
-
- # create diagonal vector of size (Nch,Nch)
- dig = np.diag(scale_factors)
- # now all rows of column 1 are multiplied with dig(1,1)
- sig = np.dot(sig,dig)
- self.sig = sig
- # 'file name;' + 'lnr;msg;'*(len(MsgList)) + '\n'
- except:
- self.sig = []
- self.error_msg = 'error: reading binary file failed'
- print('========================================================')
- print(self.error_msg)
- print(self.file_path)
- print(self.file_name)
- print('========================================================')
-
- def read_ascii(self, usecols=None):
+ if self.debug:
+ stop = time() - start
+ print('time to load HAWC2 file:', stop, 's')
- try:
- go_ascii = os.path.join(self.file_path, self.file_name) + '.dat'
-# self.sig = np.genfromtxt(go_ascii)
- self.sig = np.loadtxt(go_ascii, usecols=usecols)
-# self.sig = np.fromfile(go_ascii, dtype=np.float32, sep=' ')
-# self.sig = self.sig.reshape((self.N, self.Nch))
- except:
- self.sig = []
- self.error_msg = 'error: reading ascii file failed'
- print('========================================================')
- print(self.error_msg)
- print(self.file_path)
- print(self.file_name)
- print('========================================================')
-
-# print '========================================================'
-# print 'ASCII reading not implemented yet'
-# print '========================================================'
-# self.sig = []
-# self.error_msg = 'error: ASCII reading not implemented yet'
def reformat_sig_details(self):
"""Change HAWC2 output description of the channels short descriptive
@@ -345,7 +158,7 @@ class LoadResults(object):
# CONFIGURATION: mappings between HAWC2 and short good output:
change_list = []
- change_list.append( ['original','new improved'] )
+ change_list.append( ['original', 'new improved'] )
# change_list.append( ['Mx coo: hub1','blade1 root bending: flap'] )
# change_list.append( ['My coo: hub1','blade1 root bending: edge'] )
@@ -359,41 +172,41 @@ class LoadResults(object):
# change_list.append( ['My coo: hub3','blade3 root bending: edge'] )
# change_list.append( ['Mz coo: hub3','blade3 root bending: torsion'] )
- change_list.append( ['Mx coo: blade1','blade1 flap'] )
- change_list.append( ['My coo: blade1','blade1 edge'] )
- change_list.append( ['Mz coo: blade1','blade1 torsion'] )
+ change_list.append(['Mx coo: blade1', 'blade1 flap'])
+ change_list.append(['My coo: blade1', 'blade1 edge'])
+ change_list.append(['Mz coo: blade1', 'blade1 torsion'])
- change_list.append( ['Mx coo: blade2','blade2 flap'] )
- change_list.append( ['My coo: blade2','blade2 edge'] )
- change_list.append( ['Mz coo: blade2','blade2 torsion'] )
+ change_list.append(['Mx coo: blade2', 'blade2 flap'])
+ change_list.append(['My coo: blade2', 'blade2 edge'])
+ change_list.append(['Mz coo: blade2', 'blade2 torsion'])
- change_list.append( ['Mx coo: blade3','blade3 flap'] )
- change_list.append( ['My coo: blade3','blade3 edeg'] )
- change_list.append( ['Mz coo: blade3','blade3 torsion'] )
+ change_list.append(['Mx coo: blade3', 'blade3 flap'])
+ change_list.append(['My coo: blade3', 'blade3 edeg'])
+ change_list.append(['Mz coo: blade3', 'blade3 torsion'])
- change_list.append( ['Mx coo: hub1','blade1 out-of-plane'] )
- change_list.append( ['My coo: hub1','blade1 in-plane'] )
- change_list.append( ['Mz coo: hub1','blade1 torsion'] )
+ change_list.append(['Mx coo: hub1', 'blade1 out-of-plane'])
+ change_list.append(['My coo: hub1', 'blade1 in-plane'])
+ change_list.append(['Mz coo: hub1', 'blade1 torsion'])
- change_list.append( ['Mx coo: hub2','blade2 out-of-plane'] )
- change_list.append( ['My coo: hub2','blade2 in-plane'] )
- change_list.append( ['Mz coo: hub2','blade2 torsion'] )
+ change_list.append(['Mx coo: hub2', 'blade2 out-of-plane'])
+ change_list.append(['My coo: hub2', 'blade2 in-plane'])
+ change_list.append(['Mz coo: hub2', 'blade2 torsion'])
- change_list.append( ['Mx coo: hub3','blade3 out-of-plane'] )
- change_list.append( ['My coo: hub3','blade3 in-plane'] )
- change_list.append( ['Mz coo: hub3','blade3 torsion'] )
+ change_list.append(['Mx coo: hub3', 'blade3 out-of-plane'])
+ change_list.append(['My coo: hub3', 'blade3 in-plane'])
+ change_list.append(['Mz coo: hub3', 'blade3 torsion'])
# this one will create a false positive for tower node nr1
- change_list.append( ['Mx coo: tower','tower top momemt FA'] )
- change_list.append( ['My coo: tower','tower top momemt SS'] )
- change_list.append( ['Mz coo: tower','yaw-moment'] )
+ change_list.append(['Mx coo: tower', 'tower top momemt FA'])
+ change_list.append(['My coo: tower', 'tower top momemt SS'])
+ change_list.append(['Mz coo: tower', 'yaw-moment'])
- change_list.append( ['Mx coo: chasis','chasis momemt FA'] )
- change_list.append( ['My coo: chasis','yaw-moment chasis'] )
- change_list.append( ['Mz coo: chasis','chasis moment SS'] )
+ change_list.append(['Mx coo: chasis', 'chasis momemt FA'])
+ change_list.append(['My coo: chasis', 'yaw-moment chasis'])
+ change_list.append(['Mz coo: chasis', 'chasis moment SS'])
- change_list.append( ['DLL inp 2: 2','tower clearance'] )
+ change_list.append(['DLL inp 2: 2', 'tower clearance'])
- self.ch_details_new = np.ndarray(shape=(self.Nch,3),dtype='<U100')
+ self.ch_details_new = np.ndarray(shape=(self.Nch, 3), dtype='<U100')
# approach: look for a specific description and change it.
# This approach is slow, but will not fail if the channel numbers change
@@ -401,10 +214,10 @@ class LoadResults(object):
for ch in range(self.Nch):
# the change_list will always be slower, so this loop will be
# inside the bigger loop of all channels
- self.ch_details_new[ch,:] = self.ch_details[ch,:]
+ self.ch_details_new[ch, :] = self.ch_details[ch, :]
for k in range(len(change_list)):
- if change_list[k][0] == self.ch_details[ch,0]:
- self.ch_details_new[ch,0] = change_list[k][1]
+ if change_list[k][0] == self.ch_details[ch, 0]:
+ self.ch_details_new[ch, 0] = change_list[k][1]
# channel description should be unique, so delete current
# entry and stop looking in the change list
del change_list[k]
@@ -483,7 +296,7 @@ class LoadResults(object):
# some channel ID's are unique, use them
ch_unique = set(['Omega', 'Ae rot. torque', 'Ae rot. power',
- 'Ae rot. thrust', 'Time', 'Azi 1'])
+ 'Ae rot. thrust', 'Time', 'Azi 1'])
ch_aero = set(['Cl', 'Cd', 'Alfa', 'Vrel', 'Tors_e', 'Alfa'])
ch_aerogrid = set(['a_grid', 'am_grid'])
@@ -492,13 +305,13 @@ class LoadResults(object):
# 'component', 'pos', 'coord', 'sensortype', 'radius',
# 'blade_nr', 'units', 'output_type', 'io_nr', 'io', 'dll',
# 'azimuth', 'flap_nr'])
- df_dict = {col:[] for col in self.cols}
+ df_dict = {col: [] for col in self.cols}
df_dict['ch_name'] = []
# scan through all channels and see which can be converted
# to sensible unified name
for ch in range(self.Nch):
- items = self.ch_details[ch,2].split(' ')
+ items = self.ch_details[ch, 2].split(' ')
# remove empty values in the list
items = misc.remove_items(items, '')
@@ -511,24 +324,24 @@ class LoadResults(object):
# -----------------------------------------------------------------
# check for all the unique channel descriptions
if self.ch_details[ch,0].strip() in ch_unique:
- tag = self.ch_details[ch,0].strip()
+ tag = self.ch_details[ch, 0].strip()
channelinfo = {}
- channelinfo['units'] = self.ch_details[ch,1]
- channelinfo['sensortag'] = self.ch_details[ch,2]
+ channelinfo['units'] = self.ch_details[ch, 1]
+ channelinfo['sensortag'] = self.ch_details[ch, 2]
channelinfo['chi'] = ch
# -----------------------------------------------------------------
# or in the long description:
# 0 1 2 3 4 5 6 and up
# MomentMz Mbdy:blade nodenr: 5 coo: blade TAG TEXT
- elif self.ch_details[ch,2].startswith('MomentM'):
+ elif self.ch_details[ch, 2].startswith('MomentM'):
coord = items[5]
bodyname = items[1].replace('Mbdy:', '')
# set nodenr to sortable way, include leading zeros
# node numbers start with 0 at the root
nodenr = '%03i' % int(items[3])
# skip the attached the component
- #sensortype = items[0][:-2]
+ # sensortype = items[0][:-2]
# or give the sensor type the same name as in HAWC2
sensortype = 'momentvec'
component = items[0][-1:len(items[0])]
@@ -540,7 +353,7 @@ class LoadResults(object):
# and tag it
pos = 'node-%s' % nodenr
- tagitems = (coord,bodyname,pos,sensortype,component)
+ tagitems = (coord, bodyname, pos, sensortype, component)
tag = '%s-%s-%s-%s-%s' % tagitems
# save all info in the dict
channelinfo = {}
@@ -551,17 +364,17 @@ class LoadResults(object):
channelinfo['component'] = component
channelinfo['chi'] = ch
channelinfo['sensortag'] = sensortag
- channelinfo['units'] = self.ch_details[ch,1]
+ channelinfo['units'] = self.ch_details[ch, 1]
# -----------------------------------------------------------------
# 0 1 2 3 4 5 6 7 and up
# Force Fx Mbdy:blade nodenr: 2 coo: blade TAG TEXT
- elif self.ch_details[ch,2].startswith('Force'):
+ elif self.ch_details[ch, 2].startswith('Force'):
coord = items[6]
bodyname = items[2].replace('Mbdy:', '')
nodenr = '%03i' % int(items[4])
# skipe the attached the component
- #sensortype = items[0]
+ # sensortype = items[0]
# or give the sensor type the same name as in HAWC2
sensortype = 'forcevec'
component = items[1][1]
@@ -572,7 +385,7 @@ class LoadResults(object):
# and tag it
pos = 'node-%s' % nodenr
- tagitems = (coord,bodyname,pos,sensortype,component)
+ tagitems = (coord, bodyname, pos, sensortype, component)
tag = '%s-%s-%s-%s-%s' % tagitems
# save all info in the dict
channelinfo = {}
@@ -583,7 +396,7 @@ class LoadResults(object):
channelinfo['component'] = component
channelinfo['chi'] = ch
channelinfo['sensortag'] = sensortag
- channelinfo['units'] = self.ch_details[ch,1]
+ channelinfo['units'] = self.ch_details[ch, 1]
# -----------------------------------------------------------------
# 0 1 2 3 4 5 6 7 8
@@ -604,7 +417,7 @@ class LoadResults(object):
# skip the attached the component
#sensortype = ''.join(items[0:2])
# or give the sensor type the same name as in HAWC2
- tmp = self.ch_details[ch,0].split(' ')
+ tmp = self.ch_details[ch, 0].split(' ')
sensortype = tmp[0]
if sensortype.startswith('State'):
sensortype += ' ' + tmp[1]
@@ -616,7 +429,7 @@ class LoadResults(object):
# and tag it
pos = 'elem-%s-zrel-%s' % (elementnr, zrel)
- tagitems = (coord,bodyname,pos,sensortype,component)
+ tagitems = (coord, bodyname, pos, sensortype, component)
tag = '%s-%s-%s-%s-%s' % tagitems
# save all info in the dict
channelinfo = {}
@@ -627,7 +440,7 @@ class LoadResults(object):
channelinfo['component'] = component
channelinfo['chi'] = ch
channelinfo['sensortag'] = sensortag
- channelinfo['units'] = self.ch_details[ch,1]
+ channelinfo['units'] = self.ch_details[ch, 1]
# -----------------------------------------------------------------
# DLL CONTROL I/O
@@ -645,17 +458,17 @@ class LoadResults(object):
# description case 3
# 0 1 2 4
# hawc_dll :echo outvec : 1
- elif self.ch_details[ch,0].startswith('DLL'):
+ elif self.ch_details[ch, 0].startswith('DLL'):
# case 3
if items[1][0] == ':echo':
# hawc_dll named case (case 3) is polluted with colons
- items = self.ch_details[ch,2].replace(':','')
+ items = self.ch_details[ch,2].replace(':', '')
items = items.split(' ')
items = misc.remove_items(items, '')
dll = items[1]
io = items[2]
io_nr = items[3]
- tag = 'DLL-%s-%s-%s' % (dll,io,io_nr)
+ tag = 'DLL-%s-%s-%s' % (dll, io, io_nr)
sensortag = ''
# case 2: no reference to dll name
elif self.ch_details[ch,2].startswith('DLL'):
@@ -671,7 +484,7 @@ class LoadResults(object):
io = items[1]
io_nr = items[2]
sensortag = ' '.join(items[3:])
- tag = 'DLL-%s-%s-%s' % (dll,io,io_nr)
+ tag = 'DLL-%s-%s-%s' % (dll, io, io_nr)
# save all info in the dict
channelinfo = {}
@@ -680,19 +493,19 @@ class LoadResults(object):
channelinfo['io_nr'] = io_nr
channelinfo['chi'] = ch
channelinfo['sensortag'] = sensortag
- channelinfo['units'] = self.ch_details[ch,1]
+ channelinfo['units'] = self.ch_details[ch, 1]
# -----------------------------------------------------------------
# BEARING OUTPUS
# bea1 angle_speed rpm shaft_nacelle angle speed
- elif self.ch_details[ch,0].startswith('bea'):
- output_type = self.ch_details[ch,0].split(' ')[1]
+ elif self.ch_details[ch, 0].startswith('bea'):
+ output_type = self.ch_details[ch, 0].split(' ')[1]
bearing_name = items[0]
- units = self.ch_details[ch,1]
+ units = self.ch_details[ch, 1]
# there is no label option for the bearing output
# and tag it
- tag = 'bearing-%s-%s-%s' % (bearing_name,output_type,units)
+ tag = 'bearing-%s-%s-%s' % (bearing_name, output_type, units)
# save all info in the dict
channelinfo = {}
channelinfo['bearing_name'] = bearing_name
@@ -704,20 +517,20 @@ class LoadResults(object):
# AERO CL, CD, CM, VREL, ALFA, LIFT, DRAG, etc
# Cl, R= 0.5 deg Cl of blade 1 at radius 0.49
# Azi 1 deg Azimuth of blade 1
- elif self.ch_details[ch,0].split(',')[0] in ch_aero:
- dscr_list = self.ch_details[ch,2].split(' ')
+ elif self.ch_details[ch, 0].split(',')[0] in ch_aero:
+ dscr_list = self.ch_details[ch, 2].split(' ')
dscr_list = misc.remove_items(dscr_list, '')
- sensortype = self.ch_details[ch,0].split(',')[0]
+ sensortype = self.ch_details[ch, 0].split(',')[0]
radius = dscr_list[-1]
# is this always valid?
- blade_nr = self.ch_details[ch,2].split('blade ')[1][0]
+ blade_nr = self.ch_details[ch, 2].split('blade ')[1][0]
# sometimes the units for aero sensors are wrong!
- units = self.ch_details[ch,1]
+ units = self.ch_details[ch, 1]
# there is no label option
# and tag it
- tag = '%s-%s-%s' % (sensortype,blade_nr,radius)
+ tag = '%s-%s-%s' % (sensortype, blade_nr, radius)
# save all info in the dict
channelinfo = {}
channelinfo['sensortype'] = sensortype
@@ -729,14 +542,14 @@ class LoadResults(object):
# -----------------------------------------------------------------
# for the induction grid over the rotor
# a_grid, azi 0.00 r 1.74
- elif self.ch_details[ch,0].split(',')[0] in ch_aerogrid:
- items = self.ch_details[ch,0].split(',')
+ elif self.ch_details[ch, 0].split(',')[0] in ch_aerogrid:
+ items = self.ch_details[ch, 0].split(',')
sensortype = items[0]
items2 = items[1].split(' ')
items2 = misc.remove_items(items2, '')
azi = items2[1]
radius = items2[3]
- units = self.ch_details[ch,1]
+ units = self.ch_details[ch, 1]
# and tag it
tag = '%s-azi-%s-r-%s' % (sensortype,azi,radius)
# save all info in the dict
@@ -756,15 +569,15 @@ class LoadResults(object):
# Induc. Vy, blco, R= 1.4 // Induced wsp Vy of blade 1 at radius 1.37, local bl coo.
# Induc. Vz, glco, R= 1.4 // Induced wsp Vz of blade 1 at radius 1.37, global coo.
# Induc. Vx, rpco, R= 8.4 // Induced wsp Vx of blade 1 at radius 8.43, RP. coo.
- elif self.ch_details[ch,0].strip()[:5] == 'Induc':
- items = self.ch_details[ch,2].split(' ')
+ elif self.ch_details[ch, 0].strip()[:5] == 'Induc':
+ items = self.ch_details[ch, 2].split(' ')
items = misc.remove_items(items, '')
blade_nr = int(items[5])
radius = float(items[8].replace(',', ''))
- items = self.ch_details[ch,0].split(',')
+ items = self.ch_details[ch, 0].split(',')
coord = items[1].strip()
component = items[0][-2:]
- units = self.ch_details[ch,1]
+ units = self.ch_details[ch, 1]
# and tag it
rpl = (coord, blade_nr, component, radius)
tag = 'induc-%s-blade-%1i-%s-r-%03.02f' % rpl
@@ -787,8 +600,8 @@ class LoadResults(object):
# -----------------------------------------------------------------
# WATER SURFACE gl. coo, at gl. coo, x,y= 0.00, 0.00
- elif self.ch_details[ch,2].startswith('Water'):
- units = self.ch_details[ch,1]
+ elif self.ch_details[ch, 2].startswith('Water'):
+ units = self.ch_details[ch, 1]
# but remove the comma
x = items[-2][:-1]
@@ -806,10 +619,10 @@ class LoadResults(object):
# -----------------------------------------------------------------
# WIND SPEED
# WSP gl. coo.,Vx
- elif self.ch_details[ch,0].startswith('WSP gl.'):
- units = self.ch_details[ch,1]
- direction = self.ch_details[ch,0].split(',')[1]
- tmp = self.ch_details[ch,2].split('pos')[1]
+ elif self.ch_details[ch, 0].startswith('WSP gl.'):
+ units = self.ch_details[ch, 1]
+ direction = self.ch_details[ch, 0].split(',')[1]
+ tmp = self.ch_details[ch, 2].split('pos')[1]
x, y, z = tmp.split(',')
x, y, z = x.strip(), y.strip(), z.strip()
@@ -821,18 +634,16 @@ class LoadResults(object):
channelinfo['pos'] = (x, y, z)
channelinfo['units'] = units
channelinfo['chi'] = ch
- channelinfo['sensortype'] = 'windspeed'
- channelinfo['component'] = direction[1:]
# WIND SPEED AT BLADE
# 0: WSP Vx, glco, R= 61.5
# 2: Wind speed Vx of blade 1 at radius 61.52, global coo.
- elif self.ch_details[ch,0].startswith('WSP V'):
- units = self.ch_details[ch,1].strip()
- direction = self.ch_details[ch,0].split(' ')[1].strip()
- blade_nr = self.ch_details[ch,2].split('blade')[1].strip()[:2]
- radius = self.ch_details[ch,2].split('radius')[1].split(',')[0]
- coord = self.ch_details[ch,2].split(',')[1].strip()
+ elif self.ch_details[ch, 0].startswith('WSP V'):
+ units = self.ch_details[ch, 1].strip()
+ direction = self.ch_details[ch, 0].split(' ')[1].strip()
+ blade_nr = self.ch_details[ch, 2].split('blade')[1].strip()[:2]
+ radius = self.ch_details[ch, 2].split('radius')[1].split(',')[0]
+ coord = self.ch_details[ch, 2].split(',')[1].strip()
radius = radius.strip()
blade_nr = blade_nr.strip()
@@ -851,11 +662,11 @@ class LoadResults(object):
# FLAP ANGLE
# 2: Flap angle for blade 3 flap number 1
- elif self.ch_details[ch,0][:7] == 'setbeta':
- units = self.ch_details[ch,1].strip()
- blade_nr = self.ch_details[ch,2].split('blade')[1].strip()
+ elif self.ch_details[ch, 0][:7] == 'setbeta':
+ units = self.ch_details[ch, 1].strip()
+ blade_nr = self.ch_details[ch, 2].split('blade')[1].strip()
blade_nr = blade_nr.split(' ')[0].strip()
- flap_nr = self.ch_details[ch,2].split(' ')[-1].strip()
+ flap_nr = self.ch_details[ch, 2].split(' ')[-1].strip()
radius = radius.strip()
blade_nr = blade_nr.strip()
@@ -917,7 +728,7 @@ class LoadResults(object):
for ch_name, channelinfo in self.ch_dict.items():
cols.update(set(channelinfo.keys()))
- df_dict = {col:[] for col in cols}
+ df_dict = {col: [] for col in cols}
df_dict['ch_name'] = []
for ch_name, channelinfo in self.ch_dict.items():
cols_ch = set(channelinfo.keys())
@@ -931,7 +742,6 @@ class LoadResults(object):
self.ch_df = pd.DataFrame(df_dict)
self.ch_df.set_index('chi', inplace=True)
-
def _data_window(self, nr_rev=None, time=None):
"""
Based on a time interval, create a proper slice object
@@ -977,7 +787,7 @@ class LoadResults(object):
i_range = int(self.Freq*time_range)
window = [0, time_range]
# in case the first datapoint is not at 0 seconds
- i_zero = int(self.sig[0,0]*self.Freq)
+ i_zero = int(self.sig[0, 0]*self.Freq)
slice_ = np.r_[i_zero:i_range+i_zero]
zoomtype = '_nrrev_' + format(nr_rev, '1.0f') + 'rev'
@@ -990,23 +800,23 @@ class LoadResults(object):
slice_ = np.r_[i_start:i_end]
window = [time[0], time[1]]
- zoomtype = '_zoom_%1.1f-%1.1fsec' % (time[0], time[1])
+ zoomtype = '_zoom_%1.1f-%1.1fsec' % (time[0], time[1])
return slice_, window, zoomtype, time_range
# TODO: general signal method, this is not HAWC2 specific, move out
- def calc_stats(self, sig, i0=0, i1=-1):
+ def calc_stats(self, sig, i0=0, i1=None):
stats = {}
# calculate the statistics values:
- stats['max'] = sig[i0:i1,:].max(axis=0)
- stats['min'] = sig[i0:i1,:].min(axis=0)
- stats['mean'] = sig[i0:i1,:].mean(axis=0)
- stats['std'] = sig[i0:i1,:].std(axis=0)
+ stats['max'] = sig[i0:i1, :].max(axis=0)
+ stats['min'] = sig[i0:i1, :].min(axis=0)
+ stats['mean'] = sig[i0:i1, :].mean(axis=0)
+ stats['std'] = sig[i0:i1, :].std(axis=0)
stats['range'] = stats['max'] - stats['min']
- stats['absmax'] = np.absolute(sig[i0:i1,:]).max(axis=0)
- stats['rms'] = np.sqrt(np.mean(sig[i0:i1,:]*sig[i0:i1,:], axis=0))
- stats['int'] = integrate.trapz(sig[i0:i1,:], x=sig[i0:i1,0], axis=0)
+ stats['absmax'] = np.absolute(sig[i0:i1, :]).max(axis=0)
+ stats['rms'] = np.sqrt(np.mean(sig[i0:i1, :]*sig[i0:i1, :], axis=0))
+ stats['int'] = integrate.trapz(sig[i0:i1, :], x=sig[i0:i1, 0], axis=0)
return stats
# TODO: general signal method, this is not HAWC2 specific, move out
@@ -1030,45 +840,7 @@ class LoadResults(object):
Damage equivalent loads for each m value.
"""
- try:
- sig_rf = rainflow_astm(signal)
- except (TypeError) as e:
- print(e)
- return []
-
- if len(sig_rf) < 1 and not sig_rf:
- return []
-
- hist_data, x, bin_avg = rfc_hist(sig_rf, no_bins)
- m = np.atleast_1d(m)
-
- eq = []
- for i in range(len(m)):
- eq.append(np.power(np.sum(0.5 * hist_data *\
- np.power(bin_avg, m[i])) / neq, 1. / m[i]))
- return eq
-
- # TODO: general signal method, this is not HAWC2 specific, move out
- def cycle_matrix(self, signal, no_bins=46, m=[3, 4, 6, 8, 10, 12]):
-
-# import fatigue_tools.fatigue as ft
-# cycles, ampl_bin_mean, ampl_bin_edges, mean_bin_mean, mean_edges \
-# = ft.cycle_matrix(signal, ampl_bins=no_bins, mean_bins=1,
-# rainflow_func=ft.rainflow_windap)
-# # in this case eq = sum( n_i*S_i^m )
-# return [np.sum(cycles * ampl_bin_mean ** _m) for _m in m]
-
- try:
- sig_rf = rainflow_astm(signal)
- except:
- return []
-
- if len(sig_rf) < 1 and not sig_rf:
- return []
-
- hist_data, x, bin_avg = rfc_hist(sig_rf, no_bins)
- m = np.atleast_1d(m)
- return [np.sum(0.5 * hist_data * bin_avg ** _m) for _m in m]
+ return eq_load(signal, no_bins=no_bins, m=m, neq=neq)[0]
def blade_deflection(self):
"""
@@ -1077,18 +849,18 @@ class LoadResults(object):
# select all the y deflection channels
db = misc.DictDB(self.ch_dict)
- db.search({'sensortype' : 'state pos', 'component' : 'z'})
+ db.search({'sensortype': 'state pos', 'component': 'z'})
# sort the keys and save the mean values to an array/list
chiz, zvals = [], []
for key in sorted(db.dict_sel.keys()):
- zvals.append(-self.sig[:,db.dict_sel[key]['chi']].mean())
+ zvals.append(-self.sig[:, db.dict_sel[key]['chi']].mean())
chiz.append(db.dict_sel[key]['chi'])
- db.search({'sensortype' : 'state pos', 'component' : 'y'})
+ db.search({'sensortype': 'state pos', 'component': 'y'})
# sort the keys and save the mean values to an array/list
chiy, yvals = [], []
for key in sorted(db.dict_sel.keys()):
- yvals.append(self.sig[:,db.dict_sel[key]['chi']].mean())
+ yvals.append(self.sig[:, db.dict_sel[key]['chi']].mean())
chiy.append(db.dict_sel[key]['chi'])
return np.array(zvals), np.array(yvals)
@@ -1113,7 +885,7 @@ class LoadResults(object):
# and save
print('saving...', end='')
- np.savetxt(fname, self.sig[:,list(map_sorting.keys())], fmt=fmt,
+ np.savetxt(fname, self.sig[:, list(map_sorting.keys())], fmt=fmt,
delimiter=delimiter, header=delimiter.join(header))
print(fname)
@@ -1170,26 +942,26 @@ def ReadEigenBody(fname, debug=False):
"""
- #Body data for body number : 3 with the name :nacelle
- #Results: fd [Hz] fn [Hz] log.decr [%]
- #Mode nr: 1: 1.45388E-21 1.74896E-03 6.28319E+02
+ # Body data for body number : 3 with the name :nacelle
+ # Results: fd [Hz] fn [Hz] log.decr [%]
+ # Mode nr: 1: 1.45388E-21 1.74896E-03 6.28319E+02
FILE = opent(fname)
lines = FILE.readlines()
FILE.close()
- df_dict = {'Fd_hz':[], 'Fn_hz':[], 'log_decr_pct':[], 'body':[]}
+ df_dict = {'Fd_hz': [], 'Fn_hz': [], 'log_decr_pct': [], 'body': []}
for i, line in enumerate(lines):
if debug: print('line nr: %5i' % i)
# identify for which body we will read the data
if line[:25] == 'Body data for body number':
body = line.split(':')[2].rstrip().lstrip()
# remove any annoying characters
- body = body.replace('\n','').replace('\r','')
+ body = body.replace('\n', '').replace('\r', '')
if debug: print('modes for body: %s' % body)
# identify mode number and read the eigenfrequencies
elif line[:8] == 'Mode nr:':
- linelist = line.replace('\n','').replace('\r','').split(':')
- #modenr = linelist[1].rstrip().lstrip()
+ linelist = line.replace('\n', '').replace('\r', '').split(':')
+ # modenr = linelist[1].rstrip().lstrip()
# text after Mode nr can be empty
try:
eigenmodes = linelist[2].rstrip().lstrip().split(' ')
@@ -1205,12 +977,12 @@ def ReadEigenBody(fname, debug=False):
for k in eigenmodes:
if len(k) > 1:
eigmod.append(k)
- #eigenmodes = eigmod
+ # eigenmodes = eigmod
else:
eigmod = eigenmodes
# remove any trailing spaces for each element
for k in range(len(eigmod)):
- eigmod[k] = float(eigmod[k])#.lstrip().rstrip()
+ eigmod[k] = float(eigmod[k]) #.lstrip().rstrip()
df_dict['body'].append(body)
df_dict['Fd_hz'].append(eigmod[0])
@@ -1260,16 +1032,16 @@ def ReadEigenStructure(file_path, file_name, debug=False, max_modes=500):
"""
- #0 Version ID : HAWC2MB 11.3
- #1 ___________________________________________________________________
- #2 Structure eigenanalysis output
- #3 ___________________________________________________________________
- #4 Time : 13:46:59
- #5 Date : 28:11.2012
- #6 ___________________________________________________________________
- #7 Results: fd [Hz] fn [Hz] log.decr [%]
- #8 Mode nr: 1: 3.58673E+00 3.58688E+00 5.81231E+00
- # Mode nr:294: 0.00000E+00 6.72419E+09 6.28319E+02
+ # 0 Version ID : HAWC2MB 11.3
+ # 1 ___________________________________________________________________
+ # 2 Structure eigenanalysis output
+ # 3 ___________________________________________________________________
+ # 4 Time : 13:46:59
+ # 5 Date : 28:11.2012
+ # 6 ___________________________________________________________________
+ # 7 Results: fd [Hz] fn [Hz] log.decr [%]
+ # 8 Mode nr: 1: 3.58673E+00 3.58688E+00 5.81231E+00
+ # Mode nr:294: 0.00000E+00 6.72419E+09 6.28319E+02
FILE = opent(os.path.join(file_path, file_name))
lines = FILE.readlines()
@@ -1280,12 +1052,12 @@ def ReadEigenStructure(file_path, file_name, debug=False, max_modes=500):
# we now the number of modes by having the number of lines
nrofmodes = len(lines) - header_lines
- modes_arr = np.ndarray((3,nrofmodes))
+ modes_arr = np.ndarray((3, nrofmodes))
for i, line in enumerate(lines):
if i > max_modes:
# cut off the unused rest
- modes_arr = modes_arr[:,:i]
+ modes_arr = modes_arr[:, :i]
break
# ignore the header
@@ -1294,9 +1066,9 @@ def ReadEigenStructure(file_path, file_name, debug=False, max_modes=500):
# split up mode nr from the rest
parts = line.split(':')
- #modenr = int(parts[1])
+ # modenr = int(parts[1])
# get fd, fn and damping, but remove all empty items on the list
- modes_arr[:,i-header_lines]=misc.remove_items(parts[2].split(' '),'')
+ modes_arr[:, i-header_lines]=misc.remove_items(parts[2].split(' '), '')
return modes_arr
@@ -1430,7 +1202,7 @@ class UserWind(object):
t1 = np.exp(-math.sqrt(z_h / h_ME))
t2 = (z - z_h) / math.sqrt(z_h * h_ME)
- t3 = ( 1.0 - (z-z_h)/(2.0*math.sqrt(z_h*h_ME)) - (z-z_h)/(4.0*z_h) )
+ t3 = (1.0 - (z-z_h)/(2.0*math.sqrt(z_h*h_ME)) - (z-z_h)/(4.0*z_h))
return a_phi * t1 * t2 * t3
@@ -1479,9 +1251,9 @@ class UserWind(object):
# assert np.allclose(np.abs(u), np.abs(u2))
# assert np.allclose(np.abs(v), np.abs(v2))
- u_full = u[:,np.newaxis] + np.zeros((3,))[np.newaxis,:]
- v_full = v[:,np.newaxis] + np.zeros((3,))[np.newaxis,:]
- w_full = np.zeros((nr_vert,nr_hor))
+ u_full = u[:, np.newaxis] + np.zeros((3,))[np.newaxis, :]
+ v_full = v[:, np.newaxis] + np.zeros((3,))[np.newaxis, :]
+ w_full = np.zeros((nr_vert, nr_hor))
return u_full, v_full, w_full, x, z
@@ -1517,10 +1289,10 @@ class UserWind(object):
w_comp = np.genfromtxt(fname, skiprows=3+2+nr_vert*2,
skip_footer=i-3-3-nr_vert*3)
v_coord = np.genfromtxt(fname, skiprows=3+3+nr_vert*3,
- skip_footer=i-3-3-nr_vert*3-3)
+ skip_footer=i-3-3-nr_vert*3-3)
w_coord = np.genfromtxt(fname, skiprows=3+3+nr_vert*3+4,
- skip_footer=i-k)
- phi_deg = np.arctan(v_comp[:,0]/u_comp[:,0])*180.0/np.pi
+ skip_footer=i-k)
+ phi_deg = np.arctan(v_comp[:, 0]/u_comp[:, 0])*180.0/np.pi
return u_comp, v_comp, w_comp, v_coord, w_coord, phi_deg
@@ -1554,10 +1326,10 @@ class UserWind(object):
np.savetxt(fid, w, fmt=fmt_uvw, delimiter=' ')
h2 = b'# v coordinates (along the horizontal, nr_hor, 0 rotor center)'
fid.write(b'%s\n' % h2)
- np.savetxt(fid, v_coord.reshape((v_coord.size,1)), fmt=fmt_coord)
+ np.savetxt(fid, v_coord.reshape((v_coord.size, 1)), fmt=fmt_coord)
h3 = b'# w coordinates (zero is at ground level, height, nr_hor)'
fid.write(b'%s\n' % h3)
- np.savetxt(fid, w_coord.reshape((w_coord.size,1)), fmt=fmt_coord)
+ np.savetxt(fid, w_coord.reshape((w_coord.size, 1)), fmt=fmt_coord)
class WindProfiles(object):
@@ -1679,9 +1451,9 @@ class Turbulence(object):
# mean velocity components at the center of the box
v1, v2 = (shape[1]/2)-1, shape[1]/2
w1, w2 = (shape[2]/2)-1, shape[2]/2
- ucent = (u[:,v1,w1] + u[:,v1,w2] + u[:,v2,w1] + u[:,v2,w2]) / 4.0
- vcent = (v[:,v1,w1] + v[:,v1,w2] + v[:,v2,w1] + v[:,v2,w2]) / 4.0
- wcent = (w[:,v1,w1] + w[:,v1,w2] + w[:,v2,w1] + w[:,v2,w2]) / 4.0
+ ucent = (u[:, v1, w1] + u[:, v1, w2] + u[:, v2, w1] + u[:, v2, w2]) / 4.0
+ vcent = (v[:, v1, w1] + v[:, v1, w2] + v[:, v2, w1] + v[:, v2, w2]) / 4.0
+ wcent = (w[:, v1, w1] + w[:, v1, w2] + w[:, v2, w1] + w[:, v2, w2]) / 4.0
# FIXME: where is this range 351:7374 coming from?? The original script
# considered a box of lenght 8192
@@ -1707,9 +1479,9 @@ class Turbulence(object):
iv = np.zeros(shape)
iw = np.zeros(shape)
- iu[:,:,:] = (u - umean)/ustd*1000.0
- iv[:,:,:] = (v - vmean)/vstd*1000.0
- iw[:,:,:] = (w - wmean)/wstd*1000.0
+ iu[:, :, :] = (u - umean)/ustd*1000.0
+ iv[:, :, :] = (v - vmean)/vstd*1000.0
+ iw[:, :, :] = (w - wmean)/wstd*1000.0
# because MATLAB and Octave do a round when casting from float to int,
# and Python does a floor, we have to round first
@@ -1741,33 +1513,33 @@ class Turbulence(object):
iu, iv, iw = self.convert2bladed(fpath, basename, shape=shape)
fid = open(fpath + basename + '.wnd', 'wb')
- fid.write(struct.pack('h', R1)) # R1
- fid.write(struct.pack('h', R2)) # R2
- fid.write(struct.pack('i', turb)) # Turb
- fid.write(struct.pack('f', 999)) # Lat
- fid.write(struct.pack('f', 999)) # rough
- fid.write(struct.pack('f', 999)) # refh
- fid.write(struct.pack('f', longti)) # LongTi
- fid.write(struct.pack('f', latti)) # LatTi
- fid.write(struct.pack('f', vertti)) # VertTi
- fid.write(struct.pack('f', dv)) # VertGpSpace
- fid.write(struct.pack('f', dw)) # LatGpSpace
- fid.write(struct.pack('f', du)) # LongGpSpace
- fid.write(struct.pack('i', shape[0]/2)) # HalfAlong
- fid.write(struct.pack('f', mean_ws)) # meanWS
- fid.write(struct.pack('f', 999.)) # VertLongComp
- fid.write(struct.pack('f', 999.)) # LatLongComp
- fid.write(struct.pack('f', 999.)) # LongLongComp
- fid.write(struct.pack('i', 999)) # Int
- fid.write(struct.pack('i', seed)) # Seed
- fid.write(struct.pack('i', shape[1])) # VertGpNum
- fid.write(struct.pack('i', shape[2])) # LatGpNum
- fid.write(struct.pack('f', 999)) # VertLatComp
- fid.write(struct.pack('f', 999)) # LatLatComp
- fid.write(struct.pack('f', 999)) # LongLatComp
- fid.write(struct.pack('f', 999)) # VertVertComp
- fid.write(struct.pack('f', 999)) # LatVertComp
- fid.write(struct.pack('f', 999)) # LongVertComp
+ fid.write(struct.pack('h', R1)) # R1
+ fid.write(struct.pack('h', R2)) # R2
+ fid.write(struct.pack('i', turb)) # Turb
+ fid.write(struct.pack('f', 999)) # Lat
+ fid.write(struct.pack('f', 999)) # rough
+ fid.write(struct.pack('f', 999)) # refh
+ fid.write(struct.pack('f', longti)) # LongTi
+ fid.write(struct.pack('f', latti)) # LatTi
+ fid.write(struct.pack('f', vertti)) # VertTi
+ fid.write(struct.pack('f', dv)) # VertGpSpace
+ fid.write(struct.pack('f', dw)) # LatGpSpace
+ fid.write(struct.pack('f', du)) # LongGpSpace
+ fid.write(struct.pack('i', shape[0]/2)) # HalfAlong
+ fid.write(struct.pack('f', mean_ws)) # meanWS
+ fid.write(struct.pack('f', 999.)) # VertLongComp
+ fid.write(struct.pack('f', 999.)) # LatLongComp
+ fid.write(struct.pack('f', 999.)) # LongLongComp
+ fid.write(struct.pack('i', 999)) # Int
+ fid.write(struct.pack('i', seed)) # Seed
+ fid.write(struct.pack('i', shape[1])) # VertGpNum
+ fid.write(struct.pack('i', shape[2])) # LatGpNum
+ fid.write(struct.pack('f', 999)) # VertLatComp
+ fid.write(struct.pack('f', 999)) # LatLatComp
+ fid.write(struct.pack('f', 999)) # LongLatComp
+ fid.write(struct.pack('f', 999)) # VertVertComp
+ fid.write(struct.pack('f', 999)) # LatVertComp
+ fid.write(struct.pack('f', 999)) # LongVertComp
# fid.flush()
# bladed2 = np.ndarray((shape[0], shape[2], shape[1], 3), dtype=np.int16)
@@ -1783,9 +1555,9 @@ class Turbulence(object):
# re-arrange array for bladed format
bladed = np.ndarray((shape[0], shape[2], shape[1], 3), dtype=np.int16)
- bladed[:,:,:,0] = iu[:,::-1,:]
- bladed[:,:,:,1] = iv[:,::-1,:]
- bladed[:,:,:,2] = iw[:,::-1,:]
+ bladed[:, :, :, 0] = iu[:, ::-1, :]
+ bladed[:, :, :, 1] = iv[:, ::-1, :]
+ bladed[:, :, :, 2] = iw[:, ::-1, :]
bladed_swap_view = bladed.swapaxes(1,2)
bladed_swap_view.tofile(fid, format='%int16')