From 4d431cc606c0151a337df48655d80db480d853e7 Mon Sep 17 00:00:00 2001
From: "Mads M. Pedersen" <mmpe@dtu.dk>
Date: Thu, 1 Sep 2016 13:36:05 +0200
Subject: [PATCH] some stuff
---
setup.cfg | 1 +
wetb/__init__.py | 6 +--
wetb/hawc2/simulation.py | 2 +-
wetb/hawc2/st_file.py | 4 +-
wetb/utils/caching.py | 7 ++-
wetb/utils/cython_compile/cython_compile.py | 2 +-
wetb/utils/geometry.py | 56 ++++++++++++++++++++-
wetb/utils/tests/test_geometry.py | 11 +++-
wetb/wind/shear.py | 26 ++++++----
wetb/wind/tests/test_Shear.py | 6 ++-
10 files changed, 99 insertions(+), 22 deletions(-)
diff --git a/setup.cfg b/setup.cfg
index aa1d709..4e37be1 100644
--- a/setup.cfg
+++ b/setup.cfg
@@ -14,6 +14,7 @@ classifiers = Development Status :: 4 - Beta,
Programming Language :: Python :: 3,
Programming Language :: Python :: 3.3,
Programming Language :: Python :: 3.4,
+ Programming Language :: Python :: 3.5,
Environment :: Console,
Intended Audience :: Education,
Intended Audience :: Science/Research,
diff --git a/wetb/__init__.py b/wetb/__init__.py
index 24f6acf..a4c36db 100644
--- a/wetb/__init__.py
+++ b/wetb/__init__.py
@@ -4,9 +4,9 @@ from __future__ import division
from __future__ import absolute_import
from future import standard_library
standard_library.install_aliases()
-import pkg_resources
-test = "TEST"
+test = "TEST"
try:
- __version__ = pkg_resources.get_distribution(__name__).version
+ import pkg_resources
+ __version__ = pkg_resources.safe_version(pkg_resources.get_distribution(__name__).version)
except:
__version__ = 'unknown'
diff --git a/wetb/hawc2/simulation.py b/wetb/hawc2/simulation.py
index 0c6f584..589fc44 100755
--- a/wetb/hawc2/simulation.py
+++ b/wetb/hawc2/simulation.py
@@ -85,7 +85,7 @@ class Simulation(object):
if not os.path.isabs(htcfilename):
htcfilename = os.path.join(modelpath, htcfilename)
self.filename = os.path.basename(htcfilename)
- self.htcFile = H2aeroHTCFile(htcfilename)
+ self.htcFile = HTCFile(htcfilename)
self.time_stop = self.htcFile.simulation.time_stop[0]
self.hawc2exe = hawc2exe
self.copy_turbulence = copy_turbulence
diff --git a/wetb/hawc2/st_file.py b/wetb/hawc2/st_file.py
index 1b22366..f3c36d4 100644
--- a/wetb/hawc2/st_file.py
+++ b/wetb/hawc2/st_file.py
@@ -22,8 +22,8 @@ class StFile(object):
- r : curved length distance from main_body node 1 [m]
- m : mass per unit length [kg/m]
- - xm : xc2-coordinate from C1/2 to mass center [m]
- - ym : yc2-coordinate from C1/2 to mass center [m]
+ - x_cg : xc2-coordinate from C1/2 to mass center [m]
+ - y_cg : yc2-coordinate from C1/2 to mass center [m]
- rix : radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m]
- riy : radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m]
- xs : xc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
diff --git a/wetb/utils/caching.py b/wetb/utils/caching.py
index 6f5df27..84a38a6 100644
--- a/wetb/utils/caching.py
+++ b/wetb/utils/caching.py
@@ -8,6 +8,7 @@ from __future__ import print_function
from __future__ import division
from __future__ import absolute_import
from future import standard_library
+import sys
standard_library.install_aliases()
import inspect
@@ -80,7 +81,11 @@ def cache_function(f):
self.cache_attr_lst.add(name)
return getattr(self, name)
- if 'reload' in inspect.getargspec(f)[0]:
+ version = sys.version_info
+ if version >= (3,3):
+ if 'reload' in inspect.signature(f).parameters.values():
+ raise AttributeError("Functions decorated with cache_function are not allowed to take a parameter called 'reload'")
+ elif 'reload' in inspect.getargspec(f)[0]:
raise AttributeError("Functions decorated with cache_function are not allowed to take a parameter called 'reload'")
return wrap
diff --git a/wetb/utils/cython_compile/cython_compile.py b/wetb/utils/cython_compile/cython_compile.py
index 2d09be1..1e17e16 100644
--- a/wetb/utils/cython_compile/cython_compile.py
+++ b/wetb/utils/cython_compile/cython_compile.py
@@ -151,7 +151,7 @@ def cython_import(import_module_name, compiler=None):
fid.close()
# compile, import compiled module and delete temporary files
- module_relname = os.path.relpath(eval(module_name).__file__, str(os.getcwd())).replace(os.path.sep, ".")[:-3]
+ module_relname = eval(module_name).__package__ + "." + module_name #"os.path.relpath(eval(module_name).__file__, str(os.getcwd())).replace(os.path.sep, ".")[:-3]
return compile_and_cleanup(import_module_name, pyx_filename, module_relname, compiler)
return eval(module_name)
diff --git a/wetb/utils/geometry.py b/wetb/utils/geometry.py
index 05c7369..e82594a 100644
--- a/wetb/utils/geometry.py
+++ b/wetb/utils/geometry.py
@@ -162,7 +162,7 @@ def rpm2rads(rpm):
return rpm * 2 * np.pi / 60
-def abvrel2xyz(alpha, beta, vrel):
+def abvrel2xyz_old(alpha, beta, vrel):
"""Convert pitot tube alpha, beta and relative velocity to local Cartesian wind speed velocities
Parameters
@@ -203,3 +203,57 @@ def abvrel2xyz(alpha, beta, vrel):
z[m] = sign_vsz[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(beta[m])) ** 2 + 1 + (np.tan(alpha[m]) / np.tan(beta[m])) ** 2))
return x, y, z
+
+def abvrel2xyz(alpha, beta, vrel):
+ """Convert pitot tube alpha, beta and relative velocity to local Cartesian wind speed velocities
+
+ x : parallel to pitot tube, direction pitot tube root to tip, i.e. normal flow gives negative x\n
+ y : component in alpha plane
+ z : component in beta plane
+
+ For typical usage where pitot tube is mounted on leading edge:\n
+ x: Opposite rotational direction\n
+ y: Direction of mean wind\n
+ z: From blade root to tip\n
+
+
+ Parameters
+ ----------
+ alpha : array_like
+ Pitot tube angle of attack [rad]. Zero for flow towards pitot tube. Positive around z-axis. I.e.
+ negative alpha (normal flow) gives positive y component
+ beta : array_like
+ Pitot tube side slip angle [rad]. Zero for flow towards pitot tube. Positive around y-axis. I.e.
+ Positive beta (normal flow due to expansion and position in front of blade) gives positive z
+ vrel : array_like
+ Pitot tube relative velocity. Positive: flow towards pitot tube
+
+ Returns
+ -------
+ x : array_like
+ Wind component away from pitot tube (positive for postive vrel and -90<beta<90)
+ y : array_like
+ Wind component in alpha plane (positive for positive alpha)
+ z : array_like
+ Wind component in beta plane (positive for negative beta)
+ """
+ alpha = np.array(alpha, dtype=np.float)
+ beta = np.array(beta, dtype=np.float)
+ vrel = np.array(vrel, dtype=np.float)
+
+ sign_vsx = ((np.abs(beta) > np.pi / 2) * 2 - 1) # -1 for |beta| < 90, +1 for |beta|>90
+ sign_vsy = -np.sign(alpha) #- for alpha > 0
+ sign_vsz = np.sign(beta) # for beta>0
+
+
+ x = sign_vsx * np.sqrt(vrel ** 2 / (1 + np.tan(alpha) ** 2 + np.tan(beta) ** 2))
+
+ m = alpha != 0
+ y = np.zeros_like(alpha)
+ y[m] = sign_vsy[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(alpha[m])) ** 2 + 1 + (np.tan(beta[m]) / np.tan(alpha[m])) ** 2))
+
+ m = beta != 0
+ z = np.zeros_like(alpha)
+ z[m] = sign_vsz[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(beta[m])) ** 2 + 1 + (np.tan(alpha[m]) / np.tan(beta[m])) ** 2))
+
+ return np.array([x, y, z]).T
diff --git a/wetb/utils/tests/test_geometry.py b/wetb/utils/tests/test_geometry.py
index 291864c..e2fbf1f 100644
--- a/wetb/utils/tests/test_geometry.py
+++ b/wetb/utils/tests/test_geometry.py
@@ -14,7 +14,7 @@ import unittest
import wetb.gtsdf
import numpy as np
from wetb.utils.geometry import rad, deg, mean_deg, sind, cosd, std_deg, xyz2uvw, \
- wsp_dir2uv, wsp_dir_tilt2uvw, tand
+ wsp_dir2uv, wsp_dir_tilt2uvw, tand, abvrel2xyz
import os
@@ -117,6 +117,15 @@ class TestGeometry(unittest.TestCase):
np.testing.assert_array_almost_equal(wsp_dir_tilt2uvw(wsp, dir, tilt, wsp_horizontal=True, dir_ref=180), np.array([x, -y, z]), 5)
np.testing.assert_array_almost_equal(xyz2uvw(*wsp_dir_tilt2uvw(wsp, dir, tilt, wsp_horizontal=True), left_handed=False), xyz2uvw(x, y, z))
+ def test_abvrel2xyz(self):
+ abvrel = np.array([(0., 0, 1), (30, 0, 1), (-30, 0, 1), (0, 30, 1), (0, -30, 1), (30, 30, 1), (30, 30, 2)])
+ abvrel[:, :2] = rad(abvrel[:, :2])
+ for (x, y, z), (a, b, vrel) in zip(abvrel2xyz(*abvrel.T), abvrel):
+ #print (deg(a), deg(b), vrel, x, y, z)
+ self.assertAlmostEqual(a, np.arctan(y / x))
+ self.assertAlmostEqual(b, np.arctan(-z / x))
+ self.assertAlmostEqual(vrel, np.sqrt(x ** 2 + y ** 2 + z ** 2))
+
if __name__ == "__main__":
#import sys;sys.argv = ['', 'Test.test_rad']
diff --git a/wetb/wind/shear.py b/wetb/wind/shear.py
index 4916a86..dd61224 100644
--- a/wetb/wind/shear.py
+++ b/wetb/wind/shear.py
@@ -19,7 +19,7 @@ def _z_u(z_u_lst):
u = np.array([np.mean(np.array([u])[:]) for _, u in z_u_lst])
return z, u
-def power_shear(alpha, z_ref, u_ref, z):
+def power_shear(alpha, z_ref, u_ref):
"""Power shear
Parameters
@@ -28,22 +28,20 @@ def power_shear(alpha, z_ref, u_ref, z):
The alpha shear parameter
z_ref : int or float
The reference height
- z : int, float or array_like
- The heights of interest
u_ref : int or float
The wind speed of the reference height
Returns
-------
- power_shear : array-like
- Wind speeds at height z
+ power_shear : function
+ Function returning for wsp at input heights: f(height) -> wsp
Example
--------
- >>> power_shear(.5, 70, [20,50,70], 9)
+ >>> power_shear(.5, 70, 9)([20,50,70])
[ 4.81070235 7.60638829 9. ]
"""
- return u_ref * (np.array(z) / z_ref) ** alpha
+ return lambda z : u_ref * (np.array(z) / z_ref) ** alpha
def fit_power_shear(z_u_lst):
@@ -72,7 +70,7 @@ def fit_power_shear(z_u_lst):
alpha, _ = np.polyfit(np.log(z / z_hub), np.log(u / u_hub), 1)
return alpha
-def fit_power_shear_ref(z_u_lst, z_ref):
+def fit_power_shear_ref(z_u_lst, z_ref, plt=None):
"""Estimate power shear parameter, alpha, from two or more specific reference heights using polynomial fit.
Parameters
@@ -84,6 +82,8 @@ def fit_power_shear_ref(z_u_lst, z_ref):
- u_z2: Wind speeds or mean wind speeds at z2
z_ref : float or int
Reference height (hub height)
+ plt : matplotlib.pyplot (or similar) or None
+ Used to plot result if not None
Returns
-------
@@ -101,7 +101,13 @@ def fit_power_shear_ref(z_u_lst, z_ref):
alpha, u_ref = x
return np.sum([(u - u_ref * (z / z_ref) ** alpha) ** 2 for z, u in z_u_lst])
z_u_lst = [(z, np.mean(u)) for z, u in z_u_lst]
- return fmin(shear_error, (.1, 10), (z_u_lst, z_ref), disp=False)
+ alpha, u_ref = fmin(shear_error, (.1, 10), (z_u_lst, z_ref), disp=False)
+ if plt:
+ z, u = list(zip(*z_u_lst))
+ plt.plot(u, z, '.')
+ z = np.linspace(min(z), max(z), 100)
+ plt.plot(power_shear(alpha, z_ref, u_ref)(z), z)
+ return alpha, u_ref
@@ -188,7 +194,7 @@ if __name__ == '__main__':
from matplotlib.pyplot import plot, show
z = np.arange(0, 211)
for alpha, c in zip([0.00001, 1, 2], ['r', 'b', 'g']):
- u = power_shear(alpha, 120, 10, z)
+ u = power_shear(alpha, 120, 10)(z)
plot(u, z, c)
plot(u.mean(), 120, c + '.')
diff --git a/wetb/wind/tests/test_Shear.py b/wetb/wind/tests/test_Shear.py
index af1f6f0..318d1c3 100644
--- a/wetb/wind/tests/test_Shear.py
+++ b/wetb/wind/tests/test_Shear.py
@@ -98,8 +98,10 @@ class TestShear(unittest.TestCase):
wsp53 = data[:, 2]
wsp21 = data[:, 4]
-
- alpha, u_ref = fit_power_shear_ref([(85, wsp85), (21, wsp21)], 87.13333)
+ import matplotlib.pyplot as plt
+ alpha, u_ref = fit_power_shear_ref([(85, wsp85), (21, wsp21)], 87.13333, plt)
+ if 0:
+ plt.show()
self.assertAlmostEqual(alpha, .5, delta=.001)
self.assertAlmostEqual(u_ref, 9, delta=.01)
--
GitLab