mirror of
https://github.com/Discngine/fpocket.git
synced 2026-06-04 20:04:22 +08:00
102 lines
2.8 KiB
Python
102 lines
2.8 KiB
Python
#!/usr/bin/python
|
|
|
|
import re,numpy as npy
|
|
import sys,os
|
|
#----------------------------------------------------------------------------------------------------------
|
|
# modify your stuff here
|
|
|
|
if len(sys.argv)>3 :
|
|
d=sys.argv[1]
|
|
if not os.path.exists(d):
|
|
sys.exit(" ERROR : the dx file provided does not exist . Breaking up.")
|
|
else :
|
|
sys.exit(" ERROR : please provide the name of the dx file, the isovalue and the outputname :\n python extractISOPdb.py path/my_dx_file.dx outputname.pdb isovalue")
|
|
|
|
|
|
inputfile=sys.argv[1]
|
|
pathOutput=sys.argv[2]
|
|
iso_value=float(sys.argv[3])
|
|
#------------------------------------------------------------------------------------------------------------
|
|
|
|
#don't touch the rest...unless you know what you do :)
|
|
|
|
|
|
|
|
f=open(inputfile,"r")
|
|
|
|
|
|
#get the axis that shows the most variation during the trajectory, this will be the leading axis
|
|
|
|
#read the header - here is an example
|
|
header=""
|
|
tmp=f.readline()
|
|
while tmp[0]!="o" :
|
|
header= header + tmp
|
|
tmp=f.readline()
|
|
#print header
|
|
|
|
#read the grid size
|
|
r=re.compile('\w+')
|
|
gsize=r.findall(tmp)
|
|
gsize=[int(gsize[-3]),int(gsize[-2]),int(gsize[-1])]
|
|
#print gsize
|
|
|
|
#read the origin of the system
|
|
line=f.readline().split()
|
|
origin=[float(line[-3]),float(line[-2]),float(line[-1])]
|
|
#print origin
|
|
|
|
#read grid space
|
|
line=f.readline().split()
|
|
deltax=[float(line[-3]),float(line[-2]),float(line[-1])]
|
|
line=f.readline().split()
|
|
deltay=[float(line[-3]),float(line[-2]),float(line[-1])]
|
|
line=f.readline().split()
|
|
deltaz=[float(line[-3]),float(line[-2]),float(line[-1])]
|
|
|
|
|
|
#pay attention here, this assumes always orthogonal normalized space, but normally it should be ok
|
|
delta=npy.array([deltax[0],deltay[1],deltaz[2]])
|
|
|
|
#read the number of data
|
|
f.readline()
|
|
r=re.compile('\d+')
|
|
n_entries=int(r.findall(f.readline())[2])
|
|
|
|
if(n_entries!=gsize[0]*gsize[1]*gsize[2]) : sys.exit("Error reading the file. The number of expected data points does not correspond to the number of labeled data points in the header.")
|
|
|
|
#create a 3D numpy array filled up with 0
|
|
|
|
|
|
#initiate xyz counter for reading the grid data
|
|
z=0
|
|
y=0
|
|
x=0
|
|
|
|
print("Reading the grid. Depending on the number of data points you have this might take a while....")
|
|
path=open(pathOutput,"w")
|
|
|
|
counter=1
|
|
for count in range(n_entries//3) :
|
|
c=f.readline().split()
|
|
if(len(c)!=3) :
|
|
print("error reading grid data")
|
|
sys.exit("exiting the program")
|
|
for i in range(3):
|
|
if (iso_value<0 and float(c[i]) < iso_value) or (iso_value > 0 and float(c[i]) > iso_value) :
|
|
path.write('ATOM %5d C PTH 1 %8.3f%8.3f%8.3f%6.2f%6.2f\n'%(counter,origin[0]+float(x)*delta[0],origin[1]+float(y)*delta[1],origin[2]+float(z)*delta[2],0.0,0.0))
|
|
counter+=1
|
|
z+=1
|
|
if z >= gsize[2]:
|
|
z=0
|
|
y+=1
|
|
if y >=gsize[1]:
|
|
y=0
|
|
x+=1
|
|
|
|
|
|
path.close()
|
|
f.close()
|
|
|
|
print("finished writing %s"%(pathOutput))
|