Changeset 2412 in lmdz_wrf for trunk

Timestamp:

Mar 24, 2019, 5:37:29 PM (6 years ago)

Author:

lfita

Message:

Adding more functions from AR6-zones script

File:

• trunk/tools/geometry_tools.py (modified) (2 diffs)

trunk/tools/geometry_tools.py

@@ -17 +17 @@
 ####### Contents:
 # deg_deci: Function to pass from degrees [deg, minute, sec] to decimal angles [rad]
+# multi_rotate_2D: Function to rotate multiple vectors by a certain angle in the plain
 # position_sphere: Function to tranform fom a point in lon, lat deg coordinates to 
 #   cartesian coordinates over an sphere
+# rotate_2D: Function to rotate a vector by a certain angle in the plain
+# rotate_line2D: Function to rotate a line given by 2 pairs of x,y coordinates by a 
+#   certain angle in the plain
+# rotate_lines2D: Function to rotate multiple lines given by mulitple pars of x,y 
+#   coordinates by a certain angle in the plain
 # surface_sphere: Function to provide an sphere as matrix of x,y,z coordinates
 # spheric_line: Function to transform a series of locations in lon, lat coordinates 
 #   to x,y,z over an 3D spaceFunction to provide coordinates of a line  on a 3D space
+
+## Shapes/objects
+# ellipse_polar: Function to determine an ellipse from its center and polar coordinates
 
 ## Plotting
@@ -95 +104 @@
 
     return coords
+
+def rotate_2D(vector, angle):
+    """ Function to rotate a vector by a certain angle in the plain
+      vector= vector to rotate [y, x]
+      angle= angle to rotate [rad]
+    >>> rotate_2D(np.array([1.,0.]), np.pi/4.)
+    [ 0.70710678 -0.70710678]
+    """
+    fname = 'rotate_2D'
+
+    rotmat = np.zeros((2,2), dtype=np.float)
+
+    rotmat[0,0] = np.cos(angle)
+    rotmat[0,1] = -np.sin(angle)
+    rotmat[1,0] = np.sin(angle)
+    rotmat[1,1] = np.cos(angle)
+
+    rotvector = np.zeros((2), dtype=np.float)
+
+    vecv = np.zeros((2), dtype=np.float)
+
+    # Unifying vector
+    modvec = vector[0]**2+vector[1]**2
+    if modvec != 0: 
+        vecv[0] = vector[1]/modvec
+        vecv[1] = vector[0]/modvec
+
+        rotvec = np.matmul(rotmat, vecv)
+        rotvec = np.where(np.abs(rotvec) < 1.e-7, 0., rotvec)
+
+        rotvector[0] = rotvec[1]*modvec
+        rotvector[1] = rotvec[0]*modvec
+
+    return rotvector
+
+def multi_rotate_2D(vectors, angle):
+    """ Function to rotate multiple vectors by a certain angle in the plain
+      line= matrix of vectors to rotate
+      angle= angle to rotate [rad]
+    >>> square = np.zeros((4,2), dtype=np.float)
+    >>> square[0,:] = [-0.5,-0.5]
+    >>> square[1,:] = [0.5,-0.5]
+    >>> square[2,:] = [0.5,0.5]
+    >>> square[3,:] = [-0.5,0.5]
+    >>> multi_rotate_2D(square, np.pi/4.)
+    [[-0.70710678  0.        ]
+     [ 0.         -0.70710678]
+     [ 0.70710678  0.        ]
+     [ 0.          0.70710678]]
+    """
+    fname = 'multi_rotate_2D'
+
+    rotvecs = np.zeros(vectors.shape, dtype=np.float)
+
+    Nvecs = vectors.shape[0]
+    for iv in range(Nvecs):
+        rotvecs[iv,:] = rotate_2D(vectors[iv,:], angle)
+
+    return rotvecs
+
+def rotate_line2D(line, angle):
+    """ Function to rotate a line given by 2 pairs of x,y coordinates by a certain 
+          angle in the plain
+      line= line to rotate as couple of points [[y0,x0], [y1,x1]]
+      angle= angle to rotate [rad]
+    >>> rotate_line2D(np.array([[0.,0.], [1.,0.]]), np.pi/4.)
+    [[ 0.          0.        ]
+     [0.70710678  -0.70710678]]
+    """
+    fname = 'rotate_2D'
+
+    rotline = np.zeros((2,2), dtype=np.float)
+    rotline[0,:] = rotate_2D(line[0,:], angle)
+    rotline[1,:] = rotate_2D(line[1,:], angle)
+
+    return rotline
+
+def rotate_lines2D(lines, angle):
+    """ Function to rotate multiple lines given by mulitple pars of x,y coordinates  
+          by a certain angle in the plain
+      line= matrix of N couples of points [N, [y0,x0], [y1,x1]]
+      angle= angle to rotate [rad]
+    >>> square = np.zeros((4,2,2), dtype=np.float)
+    >>> square[0,0,:] = [-0.5,-0.5]
+    >>> square[0,1,:] = [0.5,-0.5]
+    >>> square[1,0,:] = [0.5,-0.5]
+    >>> square[1,1,:] = [0.5,0.5]
+    >>> square[2,0,:] = [0.5,0.5]
+    >>> square[2,1,:] = [-0.5,0.5]
+    >>> square[3,0,:] = [-0.5,0.5]
+    >>> square[3,1,:] = [-0.5,-0.5]
+    >>> rotate_lines2D(square, np.pi/4.)
+    [[[-0.70710678  0.        ]
+      [ 0.         -0.70710678]]
+
+     [[ 0.         -0.70710678]
+      [ 0.70710678  0.        ]]
+
+     [[ 0.70710678  0.        ]
+      [ 0.          0.70710678]]
+
+     [[ 0.          0.70710678]
+      [-0.70710678  0.        ]]]
+    """
+    fname = 'rotate_lines2D'
+
+    rotlines = np.zeros(lines.shape, dtype=np.float)
+
+    Nlines = lines.shape[0]
+    for il in range(Nlines):
+        line = np.zeros((2,2), dtype=np.float)
+        line[0,:] = lines[il,0,:]
+        line[1,:] = lines[il,1,:]
+
+        rotlines[il,:,:] = rotate_line2D(line, angle)
+
+    return rotlines
+
+####### ###### ##### #### ### ## #
+# Shapes/objects
+
+def ellipse_polar(c, a, b, Nang=100):
+    """ Function to determine an ellipse from its center and polar coordinates
+        FROM: https://en.wikipedia.org/wiki/Ellipse
+      c= coordinates of the center
+      a= distance major axis
+      b= distance minor axis
+      Nang= number of angles to use
+    """
+    fname = 'ellipse_polar'
+
+    if np.mod(Nang,2) == 0: Nang=Nang+1
+  
+    dtheta = 2*np.pi/(Nang-1)
+
+    ellipse = np.zeros((Nang,2), dtype=np.float)
+    for ia in range(Nang):
+        theta = dtheta*ia
+        rad = a*b/np.sqrt( (b*np.cos(theta))**2 + (a*np.sin(theta))**2 )
+        x = rad*np.cos(theta)
+        y = rad*np.sin(theta)
+        ellipse[ia,:] = [y+c[0],x+c[1]]
+
+    return ellipse
+
 
 ####### ####### ##### #### ### ## #