PVOX-自定义函数readoutput分析
过冷水诚挚邀请你加入Matlab仿真秀官方交流群进行Matlab学习、问题咨询、 Matlab相关资料下载，群号：927550334
“庖丁为文惠君解牛，手之所触，肩之所倚，足之所履，膝之所踦，砉然向然，奏刀騞然，莫不中音。合于《桑林》之舞，乃中《经首》之会”。
牛虽难解，十载不费一刀，徒心中有牛。复杂工程不要慌，看过冷水带你一步一步解析PVOX工具包，本期先看看自定义函数readoutput函数的构建。
[data, wfn, D_valid, W_valid] = readOutput(outName, datName)
（1）自定义函数readOutput（）的作用是读取后缀为*.out和*.dat的两个文件，输出data、wfn、D_valid、W_valid对象；
（2）Data：的作用是提取*.out中的一些变量和对应的数据，重新储存在data对象中。
（3）wfn：是提取后缀为*.dat的文件中的数据，储存在wfn中。表征波函数
（4）D_valid、W_valid这两个量是用于监控Data、wfn过程环节是否出错而设置的。假设wfn计算没有问题则W_valid赋值为1，否则为0；
outName='parsec_grid0_4.out'
datName='parsec_grid0_4.dat'
D_valid = false;
W_valid = [0,0,0,0,0];
data = [];
wfn = [];
该步骤为导入两个文件。定义输出对象。需要注意的是:D_valid为单值，W_valid为多值,这是因为存储wfn数据过程中有多个子环节需要判断。
接下是try-catch-end结构
try
Statements1
catch
Statements2
end
该结构的含义是尝试执行命令。先执行try层下的语句命令sta1，若是正常执行，则该结构功能结束，若是try层语句命令不能正常执行，则执行catch 层下的语句命令sta2，如果sta1、sta2皆不能正常执行，则跳过该结构代码，执行后续命令，不报错。这样来配合变量监控就可以从整体来把控代码有几处问题了，而不是代码以有错，就停止判断，无法了解后续代码情况。
h = waitbar(0,'Opening parsec.out');
[str,maxsize,endian] = computer;
fName = fopen('pTemp/fname.dat','wt');
fName = fopen('pTemp/fname.dat','wt');
fwrite(fName,datName);
fclose(fName);
waitbar(0.1,h);
该部分的语句的含义是新建一个文件，将后缀为*.out和*.dat的两个文件的路径写入该文件中，以便后续使用。
原始代码的第一个错误就在这里出现了。我们以截图的形式比较两段代码的差别。
过冷水的修改代码中多了两行命令，
fprintf(fName,'\n');
fwrite(fName,datName);
What you want do!过冷水总会为自己擅自修改代码编造出合理的理由，我们结合这两段代码来看另一段代码。
如果不对matlab中fname.dat写入文件时进行适当修改，下图的代码2就会报错。代码1、2是python码，不了解python的加入我们python课程，懂一点python让我们的距离更加近一点！在此说明一下:
fname = file('./pTemp/fname.dat', 'r')
fileName = fname.readline()
目的为打开fname.dat文件，并读取第一行的内容。我们们将*.out文件路径写入fname.dat中，让其读取，没毛病！
fname = file('./pTemp/fname.dat', 'r')
fileName = fname.readline()
fileName = fname.readline()
switchEndian = fname.readline()
inFile = file(fileName, 'rb')
目的为打开fname.dat文件，并读取第二行、第三行的内容。什么时候往fname.dat文件中第二行和第三行写入数据了？能不报错吗？这就是为什么添加写入内容的原因。这里改写需要注意一下各种error！
想要给读者展示错误代码报错，容易的很，在此过冷水只展示了几种没有明显错误，但实际是错的案例。过冷水只是想往文本里多添一条绝对路径，需要注意的点就有这么多，可见代码的编写细节很多，不断学习才能够完善编程知识。需要你精通matlab的跟着过冷水一行一行看代码！
if strncmp(str,'PCWIN',5)
try % Parsec Abridged Output
[s,w] = dos('del /Q pTemp\main.dat');
!"./pTools/PARSEC_PARSER_ABRIDGED.py"
catch % Parsec Un-Abridged Output
!"./pTools/PARSEC_PARSER.py"
end
else:
[s,w] = unix('rm -f pTemp/main.dat');
新知识[!"fName"]该格式是在matlab中打开文件的命令。
该段代码是判断我们运行matlab的系统是什么，不同的系统用不同的方式打开*.py文件。继续接着走：
inMain = load('pTemp/main.dat');
data.timeDiag = inMain(:,1);
data.FermiLevel = inMain(:,2);
data.timeHart = inMain(:,3);
data.ChargeDmin = inMain(:,4);
data.ChargeDmax = inMain(:,5);
data.Eig_Energy = inMain(:,6);
data.H_Energy = inMain(:,7);
data.Int_Vxc_rho = inMain(:,8);
data.Exc = inMain(:,9);
data.E_I_Energy = inMain(:,10);
data.I_I_Energy = inMain(:,11);
data.Enew_Eold = inMain(:,12);
data.Tot_Energy = inMain(:,13);
data.Energy_atom = inMain(:,14);
data.SREpot = inMain(:,15);
data.ChrgWghtPot = inMain(:,16);
这段代码的目的为从main.dat(运行*.py文件得到)中导入数据，储存到data对象中。这里再次get一个新知识点。以前过冷水就只知道：
data.timeDiag = inMain(:,1)
data= inMain(:,1)
现在变了一种赋值形式：
>> inMain = load('pTemp/main.dat');
data.timeDiag = inMain(:,1)
data =
timeDiag: [9x1 double]
Dat.char这种赋值方式很有用的，可以让我们很容易储存不同类型的数据到一个对象中，直接表示数据类型，不需要做其它操作，类似于元包组。
fOptions = fopen('pTemp/misc.dat');
opt = fgetl(fOptions);
data.num_states = str2num(strtok(opt));
opt = fgetl(fOptions);
data.num_atoms = str2num(strtok(opt));
iter = size(data.timeDiag);
data.num_iterations = iter(1);
这段代码有三个新函数：fgetl（）、str2num（）；
fgetl（）：该函数的目的是每次运行从misc.dat文件中获取一行内容。fgetl()有记忆功能。这里要注意运行次数,稍微运行次数出错会导致写入的数据对不上号，其实这里如果可以用正则匹配或者关键字定位行就不容易出错了；
strtok()：函数的含义是从字符串中找出数值字符串；
str2num：函数的目的是将字符串转化为数值。
inEigs = load('pTemp/eigenvalues.dat');
m = data.num_states;
n = size(inEigs(:,1));
n = n(1);
n = n/m;
data.State = reshape(inEigs(:,1),m,n);
data.EigenvalueRY = reshape(inEigs(:,2),m,n);
data.EigenvalueEV = reshape(inEigs(:,3),m,n);
data.Occup = reshape(inEigs(:,4),m,n);
data.Repr = reshape(inEigs(:,5),m,n)
这段代码就没有什么新意了，打开一个文件读取文件内容，并用一定规则重构元素排列。Reshape（）、size（）都是常用命令。
inForces = load('pTemp/forces.dat');
forceIon = reshape(inForces(1,:),3,data.num_atoms);
data.forceIon = transpose(forceIon);
forceSelfC = reshape(inForces(2,:),3,data.num_atoms);
data.forceSelfC = transpose(forceSelfC);
forceCoreC = reshape(inForces(3,:),3,data.num_atoms);
data.forceCoreC = transpose(forceCoreC);
代码调试不通，困扰过冷水的第二处错误就此出现。load()加载文件、reshape重构数据结构、transpose()函数没见过？调用出错？
函数transpose()：看上去定义又长，好像很重要的样子，过冷水一查：
B=transpose(A)=A’结束！
既然语句看上去没有问题，那么会有什么问题呢？这里过冷水怎么都没有想到这个坑，逻辑错误！什么意思呢？你要加载这个文件，首先的有这个文件吧！没有文件怎么加载。如果没有这个文件那就是error呗！
解释的有点费劲。*.out生成文件一般默认不包含forceIon、 forceSelfC、forceCoreC数据。然而这里的程序包默认*.out含有相应数据，不做有无的判断！这是程序设计的漏洞，因为这个原因，*.py程序编写的也有问题。过冷水的解决办法是用NaN 填充相应数据。这里可以让程序运行正确，我们不考虑科学合理性的问题。
inForces2 = load('pTemp/forces2.dat');
forceLocal = reshape(inForces2(1,:),6,data.num_atoms);
data.forceLocal = transpose(forceLocal);
forceNonLcl = reshape(inForces2(2,:),6,data.num_atoms);
data.forceNonLcl = transpose(forceNonLcl);
forceTotal = reshape(inForces2(3,:),6,data.num_atoms);
data.forceTotal = transpose(forceTotal);
forceTotal2 = reshape(inForces2(4,:),6,data.num_atoms);
data.forceTotal2 = transpose(forceTotal2);
inForces3 = load('pTemp/forces3.dat');
data.dipole = reshape(inForces3(1,:),1,3);
data.netForce = reshape(inForces3(2,:),1,3);
data.netTorque = reshape(inForces3(3,:),1,3);
D_valid = true;
data数据写入完毕!该段代码没有什么好讲的，和前面的功能相似，注意D_valid = true命令！这里就相当于在做是否有完整执行data数据写入的代码。该段代码让过冷水的涨见识了，简单的写数据的操作就几百行代码，可唬人了。经过冷水一番解刨，不过如此！
接下来回是wfn数据的写入，这里有涉及一个新自定义函数
[wfn, pot, rho, wavedata, XYZ] = wfnread(datName, h)
这真他*的是一个牵一发，动全身的工程。该完整函数程序比较复杂，初次学习涉及到的新函数比较多，需讲解内容较多，本期只讲data数据的提取。欲学后情，且下次再讲。
过冷水发表于 仿真秀 平台原创文章，未经授权禁止私自转载，如需转载请需要和作者沟通表明授权声明，未授权文章皆视为侵权行为，必将追责。如果您希望加入Matlab仿真秀官方交流群进行Matlab学习、问题咨询、 Matlab相关资料下载均可加群：927550334。
精品回顾
过冷水和你分享 matlab读取存储各种文件的方法 文末有独家金曲分享
