下面是从分子动力学轨迹数据中读取速度值的脚本。我有许多轨迹文件,其名称模式如下:
waters1445-MD001-run0100.traj
waters1445-MD001-run0200.traj
waters1445-MD001-run0300.traj
waters1445-MD001-run0400.traj
waters1445-MD001-run0500.traj
waters1445-MD001-run0600.traj
waters1445-MD001-run0700.traj
waters1445-MD001-run0800.traj
waters1445-MD001-run0900.traj
waters1445-MD001-run1000.traj
waters1445-MD002-run0100.traj
waters1445-MD002-run0200.traj
waters1445-MD002-run0300.traj
waters1445-MD002-run0400.traj
waters1445-MD002-run0500.traj
waters1445-MD002-run0600.traj
waters1445-MD002-run0700.traj
waters1445-MD002-run0800.traj
waters1445-MD002-run0900.traj
waters1445-MD002-run1000.traj
每个文件有 200 帧数据进行分析。所以我计划以这样一种方式,该代码应该一个接一个地读取每个 traj 文件(如上所示),并提取速度值并写入特定文件(text_file = open("Output.traj.dat", "a") 对应于各自的输入轨迹文件。
所以我定义了一个名为'loops(mmm)'的函数,其中'mmm'是函数'loops'的轨迹文件名解析器。
#!/usr/bin/env python
'''
always put #!/usr/bin/env python at the shebang
'''
#from __future__ import print_function
from Scientific.IO.NetCDF import NetCDFFile as Dataset
import itertools as itx
import sys
#####################
def loops(mmm):
inputfile = mmm
for FRAMES in range(0,200):
frame = FRAMES
text_file = open("Output.mmm.dat", "a")
def grouper(n, iterable, fillvalue=None):
args = [iter(iterable)] * n
return itx.izip_longest(fillvalue=fillvalue, *args)
formatxyz = "%12.7f%12.7f%12.7f%12.7f%12.7f%12.7f"
formatxyz_size = 6
formatxyzshort = "%12.7f%12.7f%12.7f"
formatxyzshort_size = 3
#ncfile = Dataset(inputfile, 'r')
ncfile = Dataset(ppp, 'r')
variableNames = ncfile.variables.keys()
#print variableNames
shape = ncfile.variables['coordinates'].shape
'''
do the header
'''
print 'title ' + str(frame)
text_file.write('title ' + str(frame) + '\n')
print "%5i%15.7e" % (shape[1],ncfile.variables['time'][frame])
text_file.write("%5i%15.7e" % (shape[1],ncfile.variables['time']\
[frame]) + '\n')
'''
do the velocities
'''
try:
xyz = ncfile.variables['velocities'][frame]
temp = grouper(2, xyz, "")
for i in temp:
z = tuple(itx.chain(*i))
if (len(z) == formatxyz_size):
print formatxyz % z
text_file.write(formatxyz % z + '\n')
elif (len(z) == formatxyzshort_size):
print formatxyzshort % z
text_file.write(formatxyzshort % z + '\n' )
except(KeyError):
xyz = [0] * shape[2]
xyz = [xyz] * shape[1]
temp = grouper(2, xyz, "")
for i in temp:
z = tuple(itx.chain(*i))
if (len(z) == formatxyz_size):
print formatxyz % z
elif (len(z) == formatxyzshort_size):
print formatxyzshort % z
x = ncfile.variables['cell_angles'][frame]
y = ncfile.variables['cell_lengths'][frame]
#text_file.close()
# program starts - generation of file name
for md in range(1,3):
if md < 10:
for pico in range(100,1100, 100):
if pico >= 1000:
kkk = "waters1445-MD00{0}-run{1}.traj".format(md,pico)
loops(kkk)
elif pico < 1000:
kkk = "waters1445-MD00{0}-run0{1}.traj".format(md,pico)
loops(kkk)
#print kkk
在 (# program starts - generation of file name) 行,代码应该生成文件名并相应地调用函数并提取速度并将值转储到 (text_file = open("Output.mmm.dat", "一个”)
执行此代码时,程序正在运行,但遗憾的是无法根据输入轨迹文件名生成输出文件。
我希望输出文件名是:
velo-waters1445-MD001-run0100.dat
velo-waters1445-MD001-run0200.dat
velo-waters1445-MD001-run0300.dat
velo-waters1445-MD001-run0400.dat
velo-waters1445-MD001-run0500.dat
.
.
.
我无法追踪我需要在哪里进行更改。