Changeset 2414 in lmdz_wrf for trunk/tools

Timestamp:

Mar 24, 2019, 8:35:53 PM (6 years ago)

Author:

lfita

Message:

Woring on 'ang_circ'

File:

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

trunk/tools/geometry_tools.py

@@ -205 +205 @@
     return rotlines
 
+def dist_points(ptA, ptB):
+    """ Function to provide the distance between two points
+      ptA: coordinates of the point A [yA, xA]
+      ptB: coordinates of the point B [yB, xB]
+    >>> dist_points([1.,1.], [-1.,-1.])
+    2.82842712475
+    """
+    fname = 'dist_points'
+
+    dist = np.sqrt( (ptA[0]-ptB[0])**2 + (ptA[1]-ptB[1])**2)
+
+    return dist
+
 ####### ###### ##### #### ### ## #
 # Shapes/objects
@@ -283 +296 @@
     return hyperbola_p, hyperbola_n
 
-def dist_points(ptA, ptB):
-    """ Function to provide the distance between two points
-      ptA: coordinates of the point A [yA, xA]
-      ptB: coordinates of the point B [yB, xB]
-    >>> dist_points([1.,1.], [-1.,-1.])
-    2.82842712475
-    """
-    fname = 'dist_points'
-
-    dist = np.sqrt( (ptA[0]-ptB[0])**2 + (ptA[1]-ptB[1])**2)
-
-    return dist
-
 def circ_sec(ptA, ptB, radii, Nang=100):
     """ Function union of point A and B by a section of a circle
@@ -314 +314 @@
         quit(-1)
 
-    if ptA[0] > ptB[0]:
-        tottheta = np.atan2(ptA[0]-ptB[0],ptA[1]-ptB[1])
-        yc=ptA[0]
-        if ptA[1] < ptB[0]: 
-            sign = -1
-            xc=ptA[1]+radii
-        else: 
-            sign = 1
-            xc=ptA[1]-radii
-    else:
-        tottheta = np.atan2(ptB[0]-ptA[0],ptB[1]-ptA[1])
-        yc=ptA[0]
-        if ptA[1] < ptB[0]: 
-            sign = 1
-            xc=ptA[1]+radii
-        else: 
-            sign = -1
-            xc=ptA[1]-radii
+    # Coordinate increments
+    dAB = np.abs(ptA-ptB)
+
+    # angle of the circular section joining points
+    tottheta = 2.*np.arcsin(distAB/2./radii)
+
+    # center along coincident bisection of the union
+    xcc = -radii
+    ycc = 0.
+
+    # Angle of the points
+    pttheta = np.arctan2(ptB[0]-ptA[0],ptB[1]-ptA[1])
+
+    # rotating the position of the center
+    rotang = np.pi/2.-pttheta
+
+    newcc = rotate_line2D(np.array([[ycc,xcc], [0.,0.]]), rotang)
+    print newcc
+
+    quit()
 
     dtheta = np.abs(tottheta)/(Nang-1)
     if sign == 1:
         dtheta = dtheta*(-1.)
+
+    print 'Lluis tottheta:', tottheta*180./np.pi, 'dtheta:', dtheta*180./np.pi, 'c:', xc,yc
 
     circ_sec = np.zeros((Nang,2), dtype=np.float)
@@ -343 +347 @@
         y = radii*np.sin(theta)
 
-        circ_sec[ia,:] = [x+xc,y+yc]
+        circ_sec[ia,:] = [y+yc,x+xc]
+        print ia, 'Lluis xy:', x,y, 'circ', circ_sec[ia,:]
     
     return circ_sec