"""
Rietan結晶構造解析に関する機能を提供するモジュール。
このモジュールは、Rietanプログラムに関連する様々な機能を提供します。
主な機能は、原子散乱因子(ASF)のパラメータの読み込みと計算、
Rietanプログラムの実行、およびRietan関連の入力ファイル(.ins)の生成です。
また、結晶学的な情報（空間群、格子定数、原子座標など）を処理し、
Rietanが利用できる形式に変換する機能も含まれます。
関連リンク:
:doc:`tkrietan_usage`
"""
import os
import numpy as np
from numpy import sin, cos, tan, arccos, arcsin, arctan, sqrt, exp, log, pi
import csv
import sys
from pprint import pprint
from tklib.tkobject import tkObject
from tklib.tkcrystal.tkcrystalobject import tkCrystalObject
from tklib.tkcrystal.tkspacegroup import tkSpaceGroup
from tklib.tkprogvars import ProgramDir, DBDir, AtomDBDir, RietanDir
from tklib.tkfile import tkFile
from tklib.tkinifile import tkIniFile
from tklib.tksci.tksci import a0, NA, torad, todeg, Round, Factorize, IsFactor
from tklib.tksci.tksci import degcos, degsin, degtan, acos,asin, atan, degacos, degasin, degatan
from tklib.tkutils import SplitFilePath, GetList
from tklib.tkre import DelQuote, DelSpace
import tklib.tkre
RietanProgramPath = os.path.join(RietanDir, "RIETAN_VENUS")
RietanPath = os.path.join(RietanProgramPath, "RIETAN64.exe")
SPGDBPath = os.path.join(RietanProgramPath, "SPGRA")
SPGRDBPath = os.path.join(RietanProgramPath, "spgr.daf")
ASFDCPath = os.path.join(RietanProgramPath, "asfdc")
[ドキュメント]
class tkRietan(tkCrystalObject):
"""
Rietan結晶構造解析に関する機能を提供するクラス。
`tkCrystalObject` を継承し、原子散乱因子や空間群の処理、
Rietan関連ファイルの読み書き、プログラム実行機能などを拡張します。
"""
def __init__(self, **args):
"""
tkRietanクラスのコンストラクタ。
空間群データベースパスとASFDCパスを初期設定し、親クラスの更新メソッドを呼び出します。
:param args: キーワード引数。親クラス `tkObject` の初期化に使用されます。
:type args: dict
"""
self.SPGDBPath = SPGDBPath
self.ASFDCPath = ASFDCPath
self.update(**args)
def __del__(self):
"""
tkRietanクラスのデストラクタ。
現在は何も処理を行いません。
"""
pass
def __str__(self):
"""
オブジェクトの文字列表現を返す。
クラスのパスを文字列として返します。
:returns: クラスのパス。
:rtype: str
"""
return self.ClassPath()
[ドキュメント]
def ReadASFParameters(self, AtomName, dbpath = ASFDCPath, XraySource = 'Cu', StopByError = 1):
"""
原子散乱因子(ASF)のパラメータをデータベースから読み込む。
`ASFDCPath` または指定されたデータベースファイルから、指定された原子のASFパラメータを検索し、
辞書形式で返します。X線源の種類によって異なる異常分散因子を選択します。
:param AtomName: 検索対象の原子名。
:type AtomName: str
:param dbpath: (optional) データベースファイルのパス。デフォルトは `ASFDCPath`。
:type dbpath: str
:param XraySource: (optional) X線源の種類 ('Cu', 'Cr', 'Fe', 'Co', 'Mo', 'Ag', 'Kb')。デフォルトは 'Cu'。
:type XraySource: str
:param StopByError: (optional) エラー発生時にプログラムを終了するかどうか。1 の場合は終了、0 の場合は続行。デフォルトは 1。
:type StopByError: int
:returns: 読み込まれたASFパラメータを含む辞書。エラーが発生した場合は0。
:rtype: dict or int
"""
db = tkFile(dbpath, 'r')
if not db:
print("Error in tkcrystal.tkrietan.ReadASFParameters: "
+ "Can not read from [{}].".format(dbpath))
if StopByError:
exit()
return 0
while 1:
line = db.ReadLine()
if not line:
print("Error in tkcrystal.tkrietan.ReadASFParameters: "
+ "Can not find ASF for [{}] from [{}].".format(AtomName, dbpath))
if StopByError:
exit()
return 0
list = tkre.Split(r'\s+', line)
self.dprint(1, list)
atomname = list[0]
atomnameonly, nametype, sitetype, charge = self.SplitAtomName(atomname)
mass = list[1]
Z = list[2]
line1 = db.ReadLine()
line2 = db.ReadLine()
line3 = db.ReadLine()
if line3 is None:
print("Error in tkcrystal.tkrietan.ReadASFParameters: Can not find ASF for [{}].".format(AtomName))
if StopByError:
exit()
return 0
if atomname.upper() == AtomName.upper():
ReadAtomName = line
break
db.Close()
inf = {'atomname': atomname, 'atomnameonly': atomnameonly,
'AtomicNumber': Z, 'Mass': mass, 'charge': charge}
a1 = tkre.Split(r"\s+", line1)
a2 = tkre.Split(r"\s+", line2)
if tkre.Search(r'0\/', line2):
a3 = []
else:
a3 = tkre.Split(r"\s+", line3)
if len(a2) == 11:
a1 = [a2[0], a2[6], a2[1], a2[7], a2[2], a2[8], a2[3], a2[9], a2[4], a2[10], a2[5]]
line2 = db.ReadLine()
else:
a1 = [a1[0], a1[1], a1[2], a1[3], a1[4], a1[5], a1[6], a1[7], 0.0, 1.0, a1[8]]
A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b = GetList(a1, 12, 0.0)
Crdf1, Crdf2, Fedf1, Fedf2, Codf1, Codf2, Cudf1, Cudf2, \
Modf1, Modf2, Agdf1, Agdf2, Kbdf1, Kbdf2 = GetList(a3, 14, 0.0)
if XraySource == 'Cr' and Crdf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Crdf1, Crdf2]
elif XraySource == 'Fe' and Fedf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Fedf1, Fedf2]
elif XraySource == 'Co' and Codf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Codf1,Codf2]
elif XraySource == 'Cu' and Cudf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Cudf1,Cudf2]
elif XraySource == 'Mo' and Modf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Modf1,Modf2]
elif XraySource == 'Ag' and Agdf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Agdf1,Agdf2]
elif XraySource == 'Kb' and Kbdf2 is not None:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b, Kbdf1,Kbdf2]
else:
a1 = [A1, B1, A2, B2, A3, B3, A4, B4, A5, B5, C, b]
inf['ASFDCParameters'] = a1
inf['ASFDCAtomName'] = AtomName
self.asfdc = inf
return inf
# $s = sin Q / lambda (in nm)
[ドキュメント]
def asf(self, s, StopByError = 1):
"""
原子散乱因子(ASF)を計算する。
`ReadASFParameters` で読み込まれたパラメータを使用して、与えられた `s` 値（sin(θ)/λ）に対するASFを計算します。
異常分散因子が存在する場合は複素数として返します。
:param s: sin(θ)/λ の値 (単位: nm⁻¹)。
:type s: float
:param StopByError: (optional) エラー発生時にプログラムを終了するかどうか。1 の場合は終了、0 の場合は続行。デフォルトは 1。
:type StopByError: int
:returns: 計算された原子散乱因子。異常分散因子が存在する場合は複素数、そうでない場合は実数。
:rtype: complex or float
"""
s2 = 0.01 * s * s
P = self.asfdc['ASFDCParameters']
n = 11; #@$pP
asf = P[n-1]
for i in range(0, n-1, 2):
asf += P[i] * exp(-P[i+1] * s2)
if P[13] is not None:
return complex(asf + P[12], P[13])
return asf
# $s = sin Q / lambda (in nm)
[ドキュメント]
def asfElectron(self, s, StopByError = 1):
"""
電子散乱因子を計算する。
`ReadASFParameters` で読み込まれた原子の情報を基に、
与えられた `s` 値（sin(θ)/λ）に対する電子散乱因子を計算します。
:param s: sin(θ)/λ の値 (単位: nm⁻¹)。
:type s: float
:param StopByError: (optional) エラー発生時にプログラムを終了するかどうか。1 の場合は終了、0 の場合は続行。デフォルトは 1。
:type StopByError: int
:returns: 計算された電子散乱因子。
:rtype: float
"""
if s <= 0.0:
s = 1.0e-5
s2 = s * s
k = 1.0e-2 * 8.0 * pi * pi / (a0 * 1.0e9)
Z = self.asfdc['AtomicNumber']
asfXray = self.asf(s, StopByError)
return k * (Z - asfXray) / s2
# $s = sin Q / lambda (in nm)
[ドキュメント]
def HydrogenAtomicScatteringFactor(self, s):
"""
水素原子の原子散乱因子を計算する。
与えられた `s` 値（sin(θ)/λ）に対する水素原子の原子散乱因子を計算します。
:param s: sin(θ)/λ の値。
:type s: float
:returns: 計算された水素原子の原子散乱因子。
:rtype: float
"""
k = 4.0 * pi * pi * a0 * a0 * 1.0e18
a = 1.0 + k * s * s
return 1.0 / a / a
"""
sub ReadHermannMauguin
{
my ($this, $ispg, $iset) = @_;
print "SPGRDBPath: $SPGRDBPath\n";
my $in = JFile->new($SPGRDBPath, 'r') or return undef;
while(1) {
my $line = $in->ReadLine();
last if(!defined $line);
my ($iSPG, $str, $otherSPGs) = ($line =~ /\(\s*(\d+)\s*\)\s*(\S.*\S)\s+(\d+)\s*$/);
#print "line: $line\n";
#print "i: $iSPG [$str] $otherSPGs\n";
next if(!defined $iSPG);
next if($iSPG != $ispg);
my @a;
$a[0] = substr($str, 0, 12);
$a[1] = substr($str, 12, 12);
$a[2] = substr($str, 24, 12);
$a[3] = substr($str, 36, 12);
for(my $i = 0 ; $i < 4 ; $i++) {
last if(!defined $a[$i]);
Utils::DelSpace($a[$i]);
}
#print " a=", join('], [', @a), "\n";
my $HM = $a[$iset-1];
#print "Hit: $iSPG=$iSPG($iset): [$HM]\n";
return $HM;
}
}
"""
"""
if($^O eq 'linux') {
$RietanFPDir = "$HOME/rietan" ;
$RietanFPProgramPath = "$HOME/rietan/bin";
$RietanFPPath = Deps::MakePath($RietanFPProgramPath, "rietan", 0);
$TeePath = "tee";
$DBDir = $RietanFPProgramPath;
$SpaceGroupDBPath = Deps::MakePath($RietanFPProgramPath, "spgra", 0);
$SPGRDBPath = Deps::MakePath($RietanFPProgramPath, "spgr.daf", 0);
}
"""
"""
my $HOME = $ENV{'HOME'};
my $WebRootDir = "d:\\MyWebs";
my $CGIPath = Deps::MakePath($WebRootDir, "cgi-bin", 1);
Jcode::convert(\$RietanFPProgramPath, $SystemCharCode) if($SystemCharCode);
#my $RietanFPPath = Deps::MakePath($RietanFPProgramPath, "rietan20040701.exe", 0);
my $TeePath = Deps::MakePath($RietanFPProgramPath, "tee.exe", 0);
my $DBDir = ProgVars::DBDir();
my $CrystalDBDir = Deps::MakePath($WebRootDir, "Research", 0);
$CrystalDBDir = Deps::MakePath($CrystalDBDir, "CrystalStructure", 0);
my $CrystalRietanFileDir = Deps::MakePath($CrystalDBDir, "Rietan", 0);
#============================================================
# 変数等取得関数
#===============================sub Read=============================
sub ClearAll { my $this=shift;
undef $this->{'FileType'}; undef $this->{'DataArray'};
undef $this->{'DiffractionPeakArray'}; }
sub GetpDiffractionPeakArray { return shift->{'DiffractionPeakArray'}; }
sub GetnDiffractionPeakArray {
my $array = shift->GetpDiffractionPeakArray(); return scalar @$array; }
sub FileType { return shift->{'FileType'}; }
sub FileName { return shift->{'FileName'}; }
sub SetFileName { my ($this,$f)=@_; return $this->{'FileName'} = $f; }
sub DataArray { return shift->{'DataArray'}; }
sub SetDataArray
{
my ($this, $DataArray) = @_;
return $this->{'DataArray'} = $DataArray;
}
sub ReadIGORPatternFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $title = $infile->ReadLine();
Utils::DelSpace($title);
$infile->SkipTo("BEGIN");
my $IntData = new GraphData;
my @Q;
my @Int;
my $line;
my $c = 0;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my @array = split(/\s+/, $line);
my ($x, $y) = Utils::RemoveSpaceElement(@array);
#print "line: $line x=$x y=$y\n";
$Q[$c] = $x + 0.0;
$Int[$c] = $y + 0.0;
$c++;
}
$IntData->SetTitle($title);
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'Intensity'} = \@Int;
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
$infile->SkipTo("BEGIN");
# my $PeakData = new GraphData;
my @Diffraction;
$Diffraction[0] = new DiffractionPeakArray;;
my @Q2;
my @Q20;
my @d;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
# my ($x1, $x0, $y, $d) = split(/\s+/, $line);
# push(@Q2, $x1);
# push(@Q20, $x0);
# push(@d, $d);
my @a = Utils::Split("\\s+", $line);
# my $n = @a;
#print "Q=$a[0]\n";
$Diffraction[0]->Add("", $a[0]+0.0, '', 0.0);
}
$infile->SkipTo("BEGIN");
for(my $i = 0 ; $i < $Diffraction[0]->nPeaks() ; $i++) {
$line = $infile->ReadLine();
last if($line =~ /^END/);
my ($h, $k, $l) = Utils::Split("\\s+", $line);
my $s = "$h $k $l";
$Diffraction[0]->Sethkl($i, $s);
}
$infile->Close();
$this->{DiffractionPeakArray} = \@Diffraction;
#print "DP: $this,$this->{DiffractionPeakArray}\n";
# $PeakData->SetTitle($title);
# $PeakData->{'2Q'} = \@Q2;
# $PeakData->{'x0'} = $PeakData->{'2Q0'} = \@Q20;
# $PeakData->{'y0'} = $PeakData->{'d'} = \@d;
# $PeakData->CalMinMax();
# $this->{'DataArray'}->AddGraphData($PeakData);
return $filename;
}
sub ReadRietPlotFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $title = $infile->ReadLine();
Utils::DelSpace($title);
Utils::DelQuote($title);
my $line = $infile->ReadLine();
my ($nData, $Start, $Step, $nPeak) = ($line =~ /^\s*(\d+)\s+([\d\.]+)\s+([\d\.]+)\s+(\d+)/);
#print "nData=$nData 2Q=$Start, $Step\n";
my $IntData = new GraphData;
my @Q;
my @Int;
my $c = 0;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my @array = split(/\s+/, $line);
@array = Utils::RemoveSpaceElement(@array);
for(my $i = 0 ; $i < @array ; $i++) {
my $Q = $Start + $c*$Step;
$Q[$c] = $Q;
$Int[$c] = $array[$i] + 0.0;
$c++;
}
last if($c >= $nData);
}
$IntData->SetTitle($title);
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'Intensity'} = \@Int;
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
# my $PeakData = new GraphData;
my @Diffraction;
$Diffraction[0] = new DiffractionPeakArray;
my @Q2;
# $c = 0;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
# my @array = split(/\s+/, $line);
# @array = Utils::RemoveSpaceElement(@array);
my @a = Utils::Split("\\s+", $line);
for(my $i = 0 ; $i < @a ; $i++) {
# $Q2[$c] = $array[$i] + 0.0;
# $c++;
$Diffraction[0]->Add("", $a[$i]+0.0, '', 0.0);
#print "Q=$a[$i]\n";
}
}
$infile->Close();
$this->{'DiffractionPeakArray'} = \@Diffraction;
# $PeakData->SetTitle($title);
# $PeakData->{'x0'} = $PeakData->{'2Q'} = \@Q2;
# $PeakData->CalMinMax();
# $this->{'DataArray'}->AddGraphData($PeakData);
return $filename;
}
sub ReadRietanFPRietPlotFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $title = $infile->ReadLine();
Utils::DelSpace($title);
Utils::DelQuote($title);
my $line = $infile->ReadLine();
#print "l1: $line";
my ($idx1, $Step, $nPeak, $idx2) = Utils::Split("\\s+", $line);
$line = $infile->ReadLine();
#print "l2: $line";
my ($nData, $Start) = Utils::Split("\\s+", $line);
print "nData=$nData 2Q=$Start, $Step nPeak=$nPeak\n";
my $IntData = new GraphData;
my @Q;
my @Obs;
my @Cal;
my @BG;
my @Diff;
my $c = 0;
while(1) {
$line = $infile->ReadLine();
my @array = Utils::Split("\\s+", $line);
for(my $i = 0 ; $i < @array ; $i += 3) {
my $Q = $Start + $c*$Step;
$Q[$c] = $Q;
$Obs[$c] = $array[$i ] + 0.0;
$Cal[$c] = $array[$i+1] + 0.0;
$BG[$c] = $array[$i+2] + 0.0;
$Diff[$c] = $Obs[$c] - $Cal[$c];
$c++;
}
last if($c >= $nData);
}
$IntData->SetTitle($title);
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'ObsIntensity'} = \@Obs;
$IntData->{'y1'} = $IntData->{'CalIntensity'} = \@Cal;
$IntData->{'y2'} = $IntData->{'BGIntensity'} = \@BG;
$IntData->{'y3'} = $IntData->{'DiffIntensity'} = \@Diff;
$IntData->CalMinMax();
$this->{DataArray}->AddGraphData($IntData);
delete $this->{DiffractionPeakArray};
my @Diffraction;
$Diffraction[0] = new DiffractionPeakArray;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my @a = Utils::Split("\\s+", $line);
for(my $i = 0 ; $i < @a ; $i++) {
$Diffraction[0]->Add("$a[$i]", $a[$i]+0.0, '', 0.0);
}
}
$infile->Close();
$this->{DiffractionPeakArray} = \@Diffraction;
# $PeakData->SetTitle($title);
# $PeakData->{'x0'} = $PeakData->{'2Q'} = \@Q2;
# $PeakData->CalMinMax();
# $this->{'DataArray'}->AddGraphData($PeakData);
return $filename;
}
sub ReadgnuplotFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $line;
my $IntData = new GraphData;
my @Q;
my @Int;
# my @hkl;
# my @Q2a;
# my @Q2b;
# my @F;
my @Diffraction;
$Diffraction[0] = new DiffractionPeakArray;
my $c = 0;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my @array = split(/\s+/, $line);
my ($Q, $Int, $h, $k, $l, $Q2a, $Q2b, $F)
= Utils::RemoveSpaceElement(@array);
$Q[$c] = $Q + 0.0;
$Int[$c] = $Int + 0.0;
if($F) {
# $hkl[$c] = "$h $k $l";
# $Q2a[$c] = $Q2a;
# $Q2b[$c] = $Q2b;
# $F[$c] = $F;
$Diffraction[0]->Add("$h $k $l", $Q2a+0.0, '', $F+0.0);
}
$c++;
}
$infile->Close();
$IntData->SetTitle('');
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'Intensity'} = \@Int;
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
$this->{'DiffractionPeakArray'} = \@Diffraction;
# my $PeakData = new GraphData;
# $PeakData->SetTitle('');
# $PeakData->{'x0'} = $PeakData->{'2Qa'} = \@Q2a;
# $PeakData->{'2Qb'} = \@Q2a;
# $PeakData->{'hkl'} = \@hkl;
# $PeakData->{'F'} = \@F;
# $PeakData->CalMinMax();
# $this->{'DataArray'}->AddGraphData($PeakData);
return $filename;
}
sub ReadIntFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $title = $infile->ReadLine();
Utils::DelSpace($title);
my $nData = $infile->ReadLine();
my $line;
my $IntData = new GraphData;
my @Q;
my @Int;
my $c = 0;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my @array = split(/\s+/, $line);
@array = Utils::RemoveSpaceElement(@array);
$Q[$c] = $array[0]+0.0;
$Int[$c] = $array[1]+0.0;
$c++;
}
$IntData->SetTitle($title);
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'Intensity'} = \@Int;
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
$infile->Close();
return $filename;
}
sub ReadIGORFittingPatternFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $line;
my $IntData = new GraphData;
my @Q;
my @Obs;
my @Cal;
my @BG;
my @Diff;
my $c = 0;
$infile->SkipTo("BEGIN");
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my ($Q, $obs, $cal, $diff, $bg) = Utils::Split("\\s+", $line);
$Q[$c] = $Q+0.0;
$Obs[$c] = $obs+0.0;
$Cal[$c] = $cal+0.0;
$Diff[$c] = $diff+0.0;
$BG[$c] = $bg+0.0;
$c++;
}
$IntData->SetTitle('');
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'ObsIntensity'} = \@Obs;
$IntData->{'y1'} = $IntData->{'CalIntensity'} = \@Cal;
$IntData->{'y2'} = $IntData->{'BGIntensity'} = \@BG;
$IntData->{'y3'} = $IntData->{'DiffIntensity'} = \@Diff;
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
$infile->SkipTo("BEGIN");
my @Diffraction;
$Diffraction[0] = new DiffractionPeakArray;;
my @Q2;
my @Q20;
my @d;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
my @a = Utils::Split("\\s+", $line);
$Diffraction[0]->Add("", $a[0]+0.0, $a[1]+0.0, 0.0);
}
$infile->SkipTo("BEGIN");
for(my $i = 0 ; $i < $Diffraction[0]->nPeaks() ; $i++) {
$line = $infile->ReadLine();
last if($line =~ /^END/);
my ($h, $k, $l) = Utils::Split("\\s+", $line);
my $s = "$h $k $l";
$Diffraction[0]->Sethkl($i, $s);
}
$infile->Close();
$this->{DiffractionPeakArray} = \@Diffraction;
#print "*** pPeakArray=$this->{DiffractionPeakArray}\n";
return $filename;
}
sub ReadFittingPatternFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $line;
my $IntData = new GraphData;
my @Q;
my @Obs;
my @Cal;
my @BG;
my @Diff;
my @Diffraction;
my $c = 0;
while($line = $infile->ReadLine()) {
last if($line =~ /^END/);
chomp($line);
my ($l1, $l2) = ($line =~ /^(.{46})(.*)$/);
my ($Q, $obs, $cal, $bg) = Utils::Split("\\s+", $l1);
$Q[$c] = $Q+0.0;
$BG[$c] = $bg+0.0;
$Obs[$c] = $obs+0.0;
$Cal[$c] = $cal+0.0;
$Diff[$c] = $obs - $cal;
#if($c < 20) {
# print "$c: l1: $l1\n";
# print " l2: $l2\n";
#}
my $ip = 0;
while(1) {
last if(length($l2) < 30);
$line = $l2;
($l1, $l2) = ($line =~ /^(.{36})(.*)$/);
#if($c < 20) {
# print " l3: $l2\n";
#}
my ($h, $k, $l, $Q2a, $Q2b, $F)
= Utils::Split("\\s+", $l1);
if(defined $F) {
$Diffraction[$ip] = new DiffractionPeakArray
unless(defined $Diffraction[$ip]);
$Diffraction[$ip]->Add("$h $k $l", $Q2a+0.0, $Q2b+0.0, $F+0.0);
}
$ip++;
}
$c++;
}
$infile->Close();
$IntData->SetTitle('');
$IntData->{'x0'} = $IntData->{'2Q'} = \@Q;
$IntData->{'y0'} = $IntData->{'ObsIntensity'} = \@Obs;
$IntData->{'y1'} = $IntData->{'CalIntensity'} = \@Cal;
$IntData->{'y2'} = $IntData->{'BGIntensity'} = \@BG;
$IntData->{'y3'} = $IntData->{'DiffIntensity'} = \@Diff;
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
$this->{'DiffractionPeakArray'} = \@Diffraction;
return $filename;
}
sub ReadTOPASTextFile
{
my ($this, $filename) = @_;
my $infile = new JFile;
my $ret = $infile->Open($filename, "r");
return undef unless($ret);
$this->SetDataArray(new GraphDataArray);
my $line = $infile->ReadLine();
$line = $infile->ReadLine();
my $title = ($line =~ /Title:(.*)$/);
Utils::DelSpace($title);
my $pos;
while(1) {
$line = $infile->ReadLine();
#print "l:$line";
if(!defined $line or $line !~ /^!/) {
$infile->seek($pos, 0) if(defined $pos);
last;
}
$pos = $infile->tell();
}
my $nData = $infile->ReadLine();
my $IntData = new GraphData;
my @Q;
my @Int;
my $c = 0;
while($line = $infile->ReadLine()) {
Utils::DelSpace($line);
last if($line eq '');
my @array = Utils::Split("\\s+", $line);
$Q[$c] = $array[0]+0.0;
$Int[$c] = $array[1]+0.0;
$c++;
}
$IntData->SetTitle($title);
$IntData->{x0} = $IntData->{'2Q'} = \@Q;
$IntData->{y0} = $IntData->{Intensity} = \@Int;
$IntData->CalMinMax();
$this->{DataArray}->AddGraphData($IntData);
$infile->Close();
return $filename;
}
sub ReadRINT2000ASCIIFile
{
my ($this, $filename) = @_;
my ($drive, $dir, $fname, $ext, $lastdir, $filebody)
= Deps::SplitFilePath($filename);
my $in = new JFile($filename, "r");
return undef unless($in);
my $line = $in->SkipTo("\\*SAMPLE ");
my ($Sample) = ($line =~ /=\s*(\S+)/);
$Sample = $filebody unless($Sample);
Utils::DelSpace($Sample);
#print "Sample: $Sample\n";
$in->rewind();
$line = $in->SkipTo("\\*START ");
my ($start) = ($line =~ /=\s*(\S+)/);
#print "Start: $start\n";
$in->rewind();
$line = $in->SkipTo("\\*STEP ");
my ($step) = ($line =~ /=\s*(\S+)/);
#print "Step: $step\n";
$in->rewind();
$line = $in->SkipTo("\\*COUNT ");
my ($nData) = ($line =~ /=\s*(\S+)/);
#print "nData: $nData\n";
my @Q2;
my @Int;
my $c = 0;
while(!$in->eof()) {
$line = $in->ReadLine();
last if($line =~ /\\*END/);
last if($line =~ /\\*EOF/);
my @a = Utils::Split("\\,\\s*", $line);
for(my $i = 0 ; $i < @a ; $i++) {
$Q2[$c] = $start + $c * $step;
$Int[$c] = $a[$i];
$c++;
last if($c >= $nData);
}
last if($c >= $nData);
}
$nData = $c;
#print "nData: $nData\n";
$in->Close();
$this->SetDataArray(new GraphDataArray);
my $IntData = new GraphData;
$IntData->SetTitle($Sample);
$IntData->{'x0'} = $IntData->{"2Q"} = \@Q2;
$IntData->{'y0'} = $IntData->{Intensity} = \@Int;
$IntData->{'x0_Name'} = "2Theta";
$IntData->SetXCaption("2Theta");
$IntData->{'y0_Name'} = $Sample;
$IntData->SetYCaption("Intensity");
$IntData->CalMinMax();
$this->{'DataArray'}->AddGraphData($IntData);
return $filename;
}
sub Read
{
my ($this, $filename) = @_;
$this->ClearAll();
my $FileType = $this->{'FileType'} = Rietan::CheckFileType($filename);
$this->SetFileName($filename);
#print "Read: this=$this\n";
#print "aFileType: $FileType\n";
if($FileType eq "Rietan Intensity File") {
return $this->ReadIntFile($filename);
}
elsif($FileType eq "Rietan Pattern File (IGOR)") {
return $this->ReadIGORPatternFile($filename);
}
elsif($FileType eq "RietanFP Pattern File (RietPlot)") {
return $this->ReadRietanFPRietPlotFile($filename);
}
elsif($FileType eq "Rietan Pattern File (RietPlot)") {
return $this->ReadRietPlotFile($filename);
}
elsif($FileType eq "Rietan Pattern File (gnuplot)") {
return $this->ReadgnuplotFile($filename);
}
elsif($FileType eq "Rietan Fitting Pattern File (IGOR)") {
return $this->ReadIGORFittingPatternFile($filename);
}
elsif($FileType eq "Rietan Fitting Pattern File") {
return $this->ReadFittingPatternFile($filename);
}
elsif($FileType eq "Rietan Input File") {
return 1;
}
elsif($FileType eq "Rietan Output File") {
return 1;
}
elsif($FileType eq 'TOPAS Text File') {
return $this->ReadTOPASTextFile($filename);
}
elsif($FileType eq "RINT2000 ASCIIFile") {
return $this->ReadRINT2000ASCIIFile($filename);
}
return undef
}
sub SetSampleName
{
my ($this, $name) = @_;
return $this->{'SampleName'} = $name;
}
sub SampleName
{
my ($this) = @_;
return $this->{'SampleName'};
}
sub SPGDBPath
{
my ($this) = @_;
return $SpaceGroupDBPath;
}
sub SPGRDBPath
{
my ($this) = @_;
return $SPGRDBPath;
}
sub GetASFDCAtoms
{
my ($this, $asfdcfile) = @_;
$asfdcfile = $ASFDCPath if(!defined $asfdcfile);
my @Atoms;
my $in = new JFile($asfdcfile, "r");
if(!$in) {
print "Can not read [$asfdcfile].\n";
return undef;
}
while(!$in->eof()) {
my $line1 = $in->ReadLine();
$in->ReadLine();
$in->ReadLine();
last if($in->eof());
my ($n) = Utils::Split("\\s+", $line1);
push(@Atoms, $n);
}
$in->Close();
return @Atoms;
}
sub FindNearestAtom
{
my ($this, $name, $charge) = @_;
my @list;
my $NearestName = '';
my $NearestCharge = 100;
#print "name: $name, charge: $charge<br>";
# my $ASFDCPath = Deps::MakePath($DBDir, "asfdc", 0);
print "Read ASFDC [$ASFDCPath] for [$name][$charge]\n";
unless(open(IN, $ASFDCPath)) {
print "Can not read ($ASFDCPath).<br>\n";
return -1;
}
while(1) {
my $line1 = <IN>;
last if(!defined $line1);
# my $line2 = <IN>;
# my $line3 = <IN>;
my ($aname) = ($line1 =~ /^(\S+)\s/);
#print "RietanFP::FindNearestAtom: aname1: $aname\n";
if($aname =~ /^($name(\.v|[\d\+-]*))$/i) {
#print "aname: $aname ";
my ($c, $sign) = ($aname =~ /(\d*)([\+-])/);
$c = 1 if($c eq '');
#print "RietanFP::FindNearestAtom: aname2: $aname c:$c s:$sign\n";
if($sign eq '-') {
$c = -$c;
}
if(abs($charge-$NearestCharge) >= abs($charge-$c)) {
$NearestCharge = $c;
$NearestName = $aname;
#print "RietanFP::FindNearestAtom: nearestname: $NearestName, charge: $NearestCharge\n";
}
push(@list, $aname);
}
}
close(IN);
return ($NearestName, @list);
}
sub FindIdenticalSpaceGroup
{
my ($this,$SPG) = @_;
my $iSPG = $SPG->iSPG();
my $Crystal = new Crystal();
$Crystal->SetLatticeParameters(3, 4, 5, 80, 85, 75);
$Crystal->SetCSpaceGroup($SPG);
$Crystal->AddAtomSite("H1", "H", 0.1, 0.15, 0.2, 1.0);
$Crystal->ExpandCoordinates();
my @AtomList = $Crystal->GetCExpandedAtomSiteList();
my $neq = @AtomList;
print "+++*** neq=$neq<br>\n";
for(my $i = 1 ; $i < 100 ; $i++) {
my $RSPG = $this->ReadSpaceGroup($iSPG, $i);
last if($RSPG == 0);
#my $ns = $RSPG->nSymmetryOperation();
#for(my $k = 0 ; $k < $ns ; $k++) {
# my $symop = $RSPG->SymmetryOperation($k+1);
# print "sym$k: $symop<br>\n";
#}
my $RCrystal = new Crystal();
$RCrystal->SetLatticeParameters(3, 4, 5, 80, 85, 75);
$RCrystal->SetCSpaceGroup($RSPG);
$RCrystal->AddAtomSite("H1", "H", 0.1, 0.15, 0.2, 1.0);
$RCrystal->ExpandCoordinates();
my @RAtomList = $RCrystal->GetCExpandedAtomSiteList();
my $Rneq = @RAtomList;
#print "+++*** $i: Rneq=$Rneq<br>\n";
next if($Rneq != $neq);
my $IsPassed = 1;
for(my $j = 0 ; $j < $neq ; $j++) {
my $atom0 = $AtomList[$j];
my ($x0, $y0, $z0) = $atom0->Position(1);
my $IsPassed0 = 0;
for(my $l = 0 ; $l < $neq ; $l++) {
my $atom1 = $RAtomList[$l];
my ($x1, $y1, $z1) = $atom1->Position(1);
my $dis = $Crystal->GetNearestInterAtomicDistance(
$x0, $y0, $z0, $x1, $y1, $z1);
#print "j=$j, l=$l: $dis: ($x0, $y0, $z0) - ($x1, $y1, $z1)<br>\n";
if($dis < 0.01) {
$IsPassed0 = 1;
last;
}
}
if($IsPassed0 == 0) {
$IsPassed = 0;
last;
}
}
next if($IsPassed == 0);
return $RSPG;
}
return -1;
}
sub Execute
{
my ($this, $filepath, $Actgion, $Program) = @_;
$Program = 'Rietan64' if(!defined $Program);
my $LF = "<br>\n";
if($Program eq 'RietanFP') {
$RietanFPDir = ProgVars::RietanFPDir();
$RietanFPProgramPath = Deps::MakePath($RietanFPDir, ["RIETAN_VENUS"], 0);
$RietanFPPath = Deps::MakePath($RietanFPProgramPath, "RIETAN64.exe", 0);
if($^O eq 'linux') {
$RietanFPDir = "$HOME/rietan" ;
$RietanFPProgramPath = "$HOME/rietan/bin";
$RietanFPPath = Deps::MakePath($RietanFPProgramPath, "rietan", 0);
}
}
#ファイル名を（ベース名, ディレクトリ名, 拡張子）に分解
# my @filenames = fileparse($filepath, "\.[^\.]+");
# my $WorkingDir = $filenames[1];
my ($drive, $directory, $filename, $ext, $lastdir, $filebody) = Deps::SplitFilePath($filepath);
my $WorkingDir = "$drive$directory";
# my $sn = $filenames[0];
my $sn = $filebody;
my $inspath = "$sn.ins";
my $intpath = "$sn.int";
my $bkgpath = "$sn.bkg";
my $patpath = "$sn.pat";
my $gpdpath = "$sn.gpd";
my $lstpath = "$sn.lst";
my $csvpath = "$sn.csv";
print "[$patpath] was deleted.$LF" if(unlink("$WorkingDir$patpath"));
print "[$lstpath] was deleted.$LF" if(unlink("$WorkingDir$lstpath"));
print "[$csvpath] was deleted.$LF" if(unlink("$WorkingDir$csvpath"));
$WorkingDir =~ s/[\\\/]$//;
$WorkingDir = "$WorkingDir\\" if($WorkingDir =~ /:$/);
if($WorkingDir ne '' and $WorkingDir ne '.') {
print "<b>Change working directory to ($WorkingDir)</b>$LF";
unless(-d $WorkingDir) {
print " ($WorkingDir) does not exist.$LF$LF";
return 0;
}
chdir($WorkingDir);
my $CurrentDirectory = cwd();
my $cd = $CurrentDirectory;
$cd =~ s/\\/\//g;
my $wd = $WorkingDir;
$wd =~ s/\\/\//g;
if(uc $cd ne uc $wd) {
print " Can not change to ($wd).$LF";
print " Current directory is ($cd)$LF$LF";
return 0;
}
}
print "<b>Set Rietan=$RietanFPProgramPath</b>$LF";
$ENV{"Rietan"} = $RietanFPProgramPath;
my $command = "\"$RietanFPPath\" \"$sn.ins\" \"$sn.int\" \"$sn.bkg\" \"$sn.itx\" "
."\"$sn.hkl\" \"$sn.xyz\" \"$sn.fos\" \"$sn.ffe\" \"$sn.fba\" \"$sn.ffi\" \"$sn.ffo\" "
."\"$sn.vesta\" \"$sn.plt\" \"$sn.gpd\"";
#"%RIETAN%\RIETAN64.exe" "%SAMPLE%.ins" "%SAMPLE%.int" "%SAMPLE%.bkg" "%SAMPLE%.itx"
#"%SAMPLE%.hkl" "%SAMPLE%.xyz" "%SAMPLE%.fos" "%SAMPLE%.ffe" "%SAMPLE%.fba" "%SAMPLE%.ffi" "%SAMPLE%.ffo"
#"%SAMPLE%.vesta" "%SAMPLE%.plt" "%SAMPLE%.gpd" > "%SAMPLE%.lst"
# my $command = "$RietanFPPath $inspath $intpath $bkgpath $patpath";
#"%RIETAN%\rietan.exe" "%SAMPLE%.ins" "%SAMPLE%.int" "%SAMPLE%.bkg" "%SAMPLE%.pat"
#"%SAMPLE%.hkl" "%SAMPLE%.xyz" "%SAMPLE%.mem" "%SAMPLE%.ffe" "%SAMPLE%.fba"
#"%SAMPLE%.ffi" "%SAMPLE%.ffo" "%SAMPLE%.vcs" | "%RIETAN%\tee.exe" "%SAMPLE%.lst"
print "$LF";
print "<b>Execute: [$command]...</b>$LF";
print "<pre>\n";
unless(open(OUT, ">$lstpath")) {
print "<b>Error in Rietan.pm:Execute: Can not write to [$lstpath].</b>$LF";
return;
}
unless(open(IN, "$command|")) {
close(OUT);
print "<b>Error in Rietan.pm:Execute: Can not execute $RietanFPPath.</b>$LF";
return;
}
while(<IN>) {
print OUT $_;
print $_;
}
close(IN);
close(OUT);
print "</pre>\n";
print "<b>Convert $gpdpath file to CSV file ($csvpath).</b>$LF";
unless(open(IN,"<$gpdpath")) {
print "<b>Error in Rietan.pm:Execute: Can not read [$gpdpath].</b>$LF";
return;
}
unless(open(OUT,">$csvpath")) {
print "<b>Error in Rietan.pm:Execute: Can not write to [$csvpath].</b>$LF";
return;
}
my $line = <IN>; #Header;
# my $line = <IN>; #Title;
# Utils::DelSpace($line);
# print "Title: $line\n";
# $line = <IN>;
# my ($ndata, $start, $step, $npeak) = Utils::Split("\\s+", $line);
# print "nData: $ndata, Start angle: $start, Step: $step, nPeaks: $npeak$LF";
my (@x, @y);
my $c = 0;
while(1) {
$line = <IN>;
last if(!defined $line);
my @a = Utils::Split("\\s+", $line);
last if(@a <= 1);
$x[$c] = $a[0];
$y[$c] = $a[1];
$c++;
}
print OUT "2Theta,Intensity\n";
# my $count = 0;
# while(<IN>) {
# my @data = split(/\s+/, $_);
# for(my $i = 0 ; $i < @data ; $i++) {
# next if($data[$i] eq '');
# my $angle = $start + $step * $count;
# $count++;
# print OUT "$angle,$data[$i]\n";
# last if($count >= $ndata);
# }
# last if($count >= $ndata);
# }
for(my $i = 0 ; $i < $c ; $i++) {
print OUT "$x[$i],$y[$i]\n";
}
close(OUT);
close(IN);
print "$LF";
$inspath = Deps::MakePath($WorkingDir, $inspath, 0);
$intpath = Deps::MakePath($WorkingDir, $intpath, 0);
$patpath = Deps::MakePath($WorkingDir, $gpdpath, 0);
# $patpath = Deps::MakePath($WorkingDir, $patpath, 0);
$lstpath = Deps::MakePath($WorkingDir, $lstpath, 0);
$csvpath = Deps::MakePath($WorkingDir, $csvpath, 0);
return ($inspath, $intpath, $patpath, $lstpath, $csvpath);
}
sub SaveInsFile
{
my ($this, $Crystal, $filename, $IsExpandCoordinate, $IsChooseRandomly,
$IsSimulation, $StartAngle, $EndAngle, $StepAngle, $Program) = @_;
$Program = 'RietanFP' if(!defined $Program);
if(!defined $IsSimulation or lc $IsSimulation eq 'simulation' or $IsSimulation eq '1') {
$IsSimulation = 1;
}
else {
$IsSimulation = 0;
}
my $LF = "<br>\n";
print "IsSimulation: $IsSimulation\n";
print "Program $Program\n";
print "Angle: ($StartAngle - $EndAngle), $StepAngle\n";
$StartAngle = 5.0 if($StartAngle eq '');
$EndAngle = 150.0 if($EndAngle eq '');
$StepAngle = 0.02 if($StepAngle eq '');
my $RefineTemplateDir = ($Program eq 'RietanFP')? ProgVars::RietanFPTemplateDir() : ProgVars::RietanTemplateDir();
#print "Temp: $RefineTemplateDir\n";
my $TemplatePath = Deps::MakePath($RefineTemplateDir, "Refinement-Template.ins", 0);
if($IsSimulation) {
$TemplatePath = Deps::MakePath($RefineTemplateDir, "Simulation-Template.ins", 0);
}
print "Template: $TemplatePath\n";
print "Save to [$filename]\n";
print "\n";
my $Content;
unless(open(IN,"<$TemplatePath")) {
print "Can not open template $TemplatePath.$LF$LF";
return;
}
while(<IN>) {
$Content .= $_;
}
close(IN);
#Rietanに使える原子名をasfdcから取得
my @AtomTypeList = $Crystal->GetCAtomTypeList();
my $nAtomType = @AtomTypeList;
my %RietanAtomName;
print "<b>Search atom/ion names available for Rietan.</b><br>\n";
for(my $i = 0 ; $i < $nAtomType ; $i++) {
my $atom = $AtomTypeList[$i];
my $name = $atom->AtomNameOnly();
my $charge = $atom->Charge() + 0;
my $namecharge = $name;
if($charge > 0) {
$namecharge = "$name$charge+";
}
elsif($charge < 0) {
my $c = -$charge;
$namecharge = "$name$c-";
}
print "\nRietanFP::SaveInsFile: Name: [$namecharge][$name][$charge]\n";
my ($NearestName, @list) = $this->FindNearestAtom($name, $charge);
$RietanAtomName{$name} = $NearestName;
print "++ $NearestName is chosen for $namecharge from (@list)<br>\n";
}
#一致する空間群を取得
print "<b>Search space group number/origin set.</b><br>\n";
my $SPG = $Crystal->GetCSpaceGroup();
my $SPGName = $SPG->SPGName();
my $iSPG = $SPG->iSPG();
my $iSPGSet = $SPG->iSet();
#print "SPG: [$SPGName] #$iSPG-$iSPGSet\n";
my $RSPG = $this->FindIdenticalSpaceGroup($SPG);
if($RSPG <= 0) {
print "<b>Error in Rietan.pm::SaveInsFile: Identical space group was not found.</b>$LF";
return -1;
}
$SPGName = $RSPG->SPGName();
$iSPG = $RSPG->iSPG();
my $iSet = $RSPG->iSet();
my $HM = $this->ReadHermannMauguin($iSPG, $iSet);
print "+++ Name: ($SPGName) [#$iSPG - $iSet]$LF";
#置換開始
my $SampleName = $this->SampleName();
$Content =~ s/{Title}/$SampleName/sig;
$Content =~ s/{Sample}/$SampleName/sig;
$Content =~ s/{PHNAME1}/$SampleName/sig;
print "T: $TemplatePath ($StartAngle, $EndAngle, $StepAngle)\n";
$Content =~ s/{StartAngle}/$StartAngle/sig;
$Content =~ s/{EndAngle}/$EndAngle/sig;
$Content =~ s/{StepAngle}/$StepAngle/sig;
my $Atoms = '';
for(my $i = 0 ; $i < $nAtomType ; $i++) {
my $atom = $AtomTypeList[$i];
my $atomname = $atom->AtomNameOnly();
my $name = $RietanAtomName{$atomname};
$Atoms .= " '$name'";# if($Atoms !~ /$name/);
}
$Content =~ s/{Elements}/$Atoms/sig;
$Content =~ s/{Atoms}/$Atoms/sig;
# my $iSPG = $Crystal->iSPGByOutputMode();
# my $iSet = $Crystal->iSetByOutputMode();
my $RietanSPG = "A-$iSPG";
$RietanSPG = "$RietanSPG-$iSet" if($iSet ne '');
$Content =~ s/{SPG}/$RietanSPG/sig;
$Content =~ s/{SPGName}/$RietanSPG/sig;
# $Content =~ s/{SPGName}/$SPGName/sig;
$Content =~ s/{HKLM1}/$HM/sig;
my ($a,$b,$c,$alpha,$beta,$gamma) = $Crystal->LatticeParametersByOutputMode(0);
my $CELLQ;
$CELLQ = "$a $b $c $alpha $beta $gamma 0.5 0000000";
$Content =~ s/{CELLQ}/$CELLQ/sig;
my @ExpandedAtomSiteList = $Crystal->GetCExpandedAtomSiteListByOutputMode();
my $nExpandedAtomSite = @ExpandedAtomSiteList;
my $ATOMCOORDS;
my %AtomCount;
for(my $i = 0 ; $i < $nExpandedAtomSite ; $i++) {
my $atom = $ExpandedAtomSiteList[$i];
my $atomname = $atom->AtomName();
my $atomnameonly = $atom->AtomNameOnly();
my $name = $RietanAtomName{$atomnameonly};
my ($x,$y,$z) = $atom->Position(1);
my $occupancy = $atom->Occupancy();
my $i1 = $i+1;
$AtomCount{$atomname}++;
my $c = $AtomCount{$atomname};
my $label = "$atomnameonly$c";
$ATOMCOORDS .= " $label/$name $occupancy $x $y $z 0.5000 01101\n";
}
$Content =~ s/{ATOMCOORDS}/$ATOMCOORDS/sig;
unless(open(OUT,">$filename")) {
print "Can not write to $filename.$LF$LF";
return;
}
print OUT $Content;
print OUT "\n";
close(OUT);
return 1;
}
1;
"""