[cig-commits] r19138 - seismo/3D/FAULT_SOURCE/branches/new_fault_db/CUBIT

[percygalvez at geodynamics.org]

Author: percygalvez
Date: 2011-11-01 04:57:38 -0700 (Tue, 01 Nov 2011)
New Revision: 19138

Added:
   seismo/3D/FAULT_SOURCE/branches/new_fault_db/CUBIT/functions.py
Log:
python split_nodes_routine for complex faults

Added: seismo/3D/FAULT_SOURCE/branches/new_fault_db/CUBIT/functions.py
===================================================================
--- seismo/3D/FAULT_SOURCE/branches/new_fault_db/CUBIT/functions.py	                        (rev 0)
+++ seismo/3D/FAULT_SOURCE/branches/new_fault_db/CUBIT/functions.py	2011-11-01 11:57:38 UTC (rev 19138)
@@ -0,0 +1,136 @@
+#!python
+
+def read_nodes_coord_file(filein):
+    file = open(filein,'r')
+    nnodes = file.readline()
+    nodxyzlist = file.readlines()
+    x=[]
+    y=[]
+    z=[]
+    node=[]
+    for nodxyz in nodxyzlist:
+        nodxyzsplit = nodxyz.split()
+        node.append(int(nodxyzsplit[0]))
+        x.append(float(nodxyzsplit[1]))
+        y.append(float(nodxyzsplit[2]))
+        z.append(float(nodxyzsplit[3]))
+    file.close() 
+    return node,x,y,z
+
+def read_fault_nodes(filein):
+    file = open(filein,'r')
+    nodes_side1_side2 = file.readline()
+    nside_ud = nodes_side1_side2.split()
+    nnodes_u = int(nside_ud[0])
+    nnodes_d = int(nside_ud[1])
+    nodes_u=[]
+    nodes_d=[]
+
+    for i in range(nnodes_u):
+        hnodes = file.readline()
+        hnodesplit = hnodes.split()
+        nodes_u.append(int(hnodesplit[1]))
+        nodes_u.append(int(hnodesplit[2]))
+        nodes_u.append(int(hnodesplit[3]))
+        nodes_u.append(int(hnodesplit[4]))
+   
+    for i in range(nnodes_d):
+        hnodes = file.readline()
+        hnodesplit = hnodes.split()
+        nodes_d.append(int(hnodesplit[1]))
+        nodes_d.append(int(hnodesplit[2]))
+        nodes_d.append(int(hnodesplit[3]))
+        nodes_d.append(int(hnodesplit[4]))
+
+    nodes_d  = set(nodes_d)  # save unique fault nodes side u
+    nodes_u  = set(nodes_u)  # save unique fault nodes side d
+    ncommond = nodes_u&nodes_d  # glue nodes . 
+    nodes_d  = nodes_d^ncommond # saving only split nodes side d
+    nodes_u  = nodes_u^ncommond # saving only split nodes side u
+
+    file.close()
+    return nodes_u,nodes_d
+
+def nodes_coords_fault_open(file_nodes_coord,fault_file,name_out,fside1,fside2): 
+    x=[]
+    y=[]
+    z=[]
+    node=[]
+    nodes_u=[]
+    nodes_d=[]
+    txt = ''
+    txt_fault = ''
+    output_file = open(name_out,'w')
+    fsideu = open(fside1,'w')
+    fsided = open(fside2,'w')
+    node,x,y,z = read_nodes_coord_file(file_nodes_coord)
+    output_file.write('%10i\n' % len(node))
+    node_u,node_d = read_fault_nodes(fault_file)
+    for i in range(len(node)):
+        if node[i] in node_u:
+             z[i] = z[i] + 500.0
+             print "node_u",node[i]
+             txt_fault=('%10i %20f %20f %20f\n') % (node[i],x[i],y[i],z[i])
+             fsideu.write(txt_fault)
+        if node[i] in node_d:
+             z[i] = z[i] - 500.0
+             print "node_d",node[i]
+             txt_fault=('%10i %20f %20f %20f\n') % (node[i],x[i],y[i],z[i])
+             fsided.write(txt_fault)
+        txt=('%10i %20f %20f %20f\n') % (node[i],x[i],y[i],z[i])
+        output_file.write(txt)
+
+    output_file.close()
+    fsideu.close()
+    fsided.close()
+    return
+
+file_nodes_coord = 'MESH/nodes_coords_file'
+name_out = 'nodes_coords_file_open_fault'
+fault_file ='MESH/fault_file_1.dat'
+fsideu = 'fault_sideu.dat'
+fsided = 'fault_sided.dat'
+nodes_coords_fault_open(file_nodes_coord,fault_file,name_out,fsideu,fsided)
+
+
+def m2km(name_in,name_out):
+    txt =''
+    km = 1000
+    input_file  = open(name_in,'r')
+    output_file = open(name_out,'w')
+    nnodes      = input_file.readline()
+    nodes_num   = int(nnodes)
+    output_file.write('%10i\n' % nodes_num)
+    nodes_list  = input_file.readlines()
+    nodes       = []
+    for node in nodes_list:
+        nodes   = node.split()
+        node_id = int(nodes[0])
+        x= float(nodes[1])*km
+        y= float(nodes[2])*km
+        z= float(nodes[3])*km
+        txt=('%10i %20f %20f %20f\n') % (node_id,x,y,z)
+        output_file.write(txt)
+
+    output_file.close()
+    return
+
+#name_in = 'MESH/nodes_coords_file_km'
+#name_out= 'nodes_coords_file'
+#m2km(name_in,name_out)
+
+
+#import numpy as np
+#def find_nearest(array,value):
+#    idx=(np.abs(array-value)).argmin()
+#    return array[idx]
+
+#array=np.random.random(10)
+#print(array)
+# [ 0.21069679  0.61290182  0.63425412  0.84635244  0.91599191  0.00213826
+#   0.17104965  0.56874386  0.57319379  0.28719469]
+
+#value=0.5
+
+#print(find_nearest(array,value))
+