基于理论结构因子校正衍射实验数据

casjxm

实验电子结构（或电子结构晶体学）研究中，对实验结构因子（或衍射强度-结构因子的平方）的校正对获得高质量的电子密度图是非常关键的，本贴介绍基于理论结构因子的方法。这种方法是经验的，并能显著改善电子密度图。

定义：
理论结构因子(fcal)：根据晶体结构计算出理论电子密度，并根据理论电子密度计算出理论的结构因子，理论的结构因子不含有实验误差。
实验结构因子(fobs)：实验测试的，含实验误差。

经验假定：
对于一个特定衍射点hkl，如果理论结构因子fcal与实验结构因子fobs的差在3*sigma（sigma为测量误差）以内，很可能是实验测试有问题，并将该hkl的fobs更改为fcal；如果超过fcal与fobs的差超过3*sigma，这很可能是实验测试的真实信号。

实例测试：
1. 先看一下常规的IAM模型精修后的残余密度图，可以看出键电荷，N和O原子上是的孤对电子，但电子云形状不规则。

2.  通过olex的NoSpherA2精修，获得.fcf文件，里面含有fcal，fobs，sigma（实际上是结构因子的平方，这里为表述方便）。使用.fcf文件中的fcal进行常规IAM精修，从残余密度可以明显看出键电荷，N和O原子上的孤对电子，而且电子云形状规则，O上的两对孤对电子呈对称分布。
olex的NoSpherA2精修可参考：https://www.matstr.com/forum.php?mod=viewthread&tid=642&extra=


3.  根据上述经验假定，校正fobs，并采用校正后的fobs进行IAM精修，其残余密度比没有经过校正的要好得多，电子云清晰，形状规则，但与理论结构因子精修出的O原子上孤对电子呈对称分布的情况不同，这里O原子上孤对电子明显呈不对称分布，应该提现了实验的真实信号，可能由于这两对孤对电子在晶体中的化学环境不同而呈现不同的特征。


Matlab代码：
%从test.fcf中读取数据
filename = 'test.fcf';
fileID = fopen(filename, 'r'); % 打开文件
data = textscan(fileID, '%d %d %d %f %f %f %c', 'Delimiter', '\t', 'HeaderLines', 32); % 读取浮点数
fclose(fileID); % 关闭文件
h= data{1};
k=data{2};
l=data{3};
fcal=data{4};
fobs=data{5};
sigma=data{6};
hkl=[h k l];intensity=[fcal fobs sigma];
n=size(h);

%保存理论结构因子为shelx精修的hkl格式
fileID = fopen('cal.hkl','w'); % 打开文件
for i=1:n(1)
if(intensity(i,1)>1e-6)
fprintf(fileID,'%4d%4d%4d%8.2f%8.2f\n',hkl(i,1),hkl(i,2),hkl(i,3),intensity(i,1),intensity(i,3)); % 将数据写入文件
end
end
fclose(fileID);

%保存实验结构因子为shelx精修的hkl格式
fileID = fopen('obs.hkl','w'); % 打开文件
for i=1:n(1)
if(intensity(i,2)>1e-6)
fprintf(fileID,'%4d%4d%4d%8.2f%8.2f\n',hkl(i,1),hkl(i,2),hkl(i,3),intensity(i,2),intensity(i,3)); % 将数据写入文件
end
end
fclose(fileID);

%保存基于理论结构因子校正后的实验结构因子为shelx精修的hkl格式
fileID = fopen('corrrected-obs.hkl','w'); % 打开文件
for i=1:n(1)
if((abs((intensity(i,2)-intensity(i,1))/intensity(i,3))>3)&&(intensity(i,1)>0.01))
fprintf(fileID,'%4d%4d%4d%8.2f%8.2f\n',hkl(i,1),hkl(i,2),hkl(i,3),intensity(i,2),intensity(i,3)); % 将数据写入文件
elseif (intensity(i,1)>0.01)
fprintf(fileID,'%4d%4d%4d%8.2f%8.2f\n',hkl(i,1),hkl(i,2),hkl(i,3),intensity(i,1),intensity(i,3)); % 将数据写入文件
end
end
fclose(fileID);