#!/usr/bin/perl

BEGIN {
use lib 'd:/Programs/Perl/lib';
#use lib '/home/tkamiya/bin/lib';
my $BaseDir = $ENV{'TkPerlDir'};
print "\n\nBaseDir: $BaseDir\n";
@INC = ("$BaseDir/lib", "$BaseDir/GULP", @INC);
#use lib "$BaseDir/lib";
#use lib "$BaseDir/GULP";
}

use strict;
#use warnings;

use Utils;

use MyApplication;

use Crystal::CIF;
use Crystal::Crystal;
use Crystal::BandStructure;

#===============================================
# デバッグ関係変数
#===============================================
#$PrintLevelが大きいほど、情報が詳しくなる
my $PrintLevel = 0;

#===============================================
# 大域変数
#===============================================
my $HOME = $ENV{'HOME'};
my $GULPDir         = Deps::MakePath($HOME, "gulp", 0);
my $GULPLibDir      = Deps::MakePath($GULPDir, "Libraries", 0);

#===============================================
# 文字コード関係変数
#===============================================
# sjis, euc, jis, noconv
my $PrintCharCode      = Deps::PrintCharCode();
my $OSCharCode         = Deps::OSCharCode();
my $FileSystemCharCode = Deps::FileSystemCharCode();

my $LF        = Deps::LF();
my $DirSep    = Deps::DirSep();
my $RegDirSep = Deps::RegDirSep();

#===============================================
# Applicationオブジェクト作成
#===============================================
my $App = new MyApplication;
exit if($App->Initialize() < 0);

#$App->SetLF("<br>\n");
#$App->SetPrintCharCode("sjis");
#$App->SetDebug($Debug);
$App->SetDeleteHTMLFlag(1);

#===============================================
# スクリプト大域変数
#===============================================
my $InitialDirectory = Deps::GetWorkingDirectory();

#==========================================
# コマンドラインオプション読み込み
#==========================================
$App->AddArgument("--Action",    "--Action=[MakeInput|MakeDOSCSV|ModifyeHtunerInput]", '');
$App->AddArgument("--aKProduct", "--aKProduct=val [Def: 2.0]", 2.0);
$App->AddArgument("--db",        "--db=eH parameter file [Def: 'eht_parms.dat']", 'eht_parms.dat');
$App->AddArgument("--KPOINTS",   "--KPOINTS=Band k point file [Def: '*.klist']",  '*.klist');
$App->AddArgument("--IBZKPT",    "--IBZKPT=SCF k point file [Def: 'IBZKPT']",     'IBZKPT');
$App->AddArgument("--Output",    "--Output=Output file [Def: 'bind.in']",      'bind.in');
$App->AddArgument("--DebugMode", "--DebugMode: Set DebugMode", '');
exit 1 if($App->ReadArgs(1) != 1);
my $Args = $App->Args();

#==========================================
# メイン関数スタート
#==========================================

#Utils::InitHTML("Research", $WebCharSet, "_self");

my $Debug = $Args->GetGetArg("DebugMode");
$App->SetDebug($Debug);
my $Action = $Args->GetGetArg("Action");

if($Action =~ /MakeInput/i) {
	&MakeInput();
}
elsif($Action =~ /MakeDOSCSV/i) {
	&MakeDOSCSV();
}
elsif($Action =~ /ModifyeHtunerInput/i) {
	&ModifyeHtunerInput();
}
else {
	$App->print("Error: Invald Action: $Action\n");
}

#Utils::EndHTML();

exit;

#===============================================
# スクリプト終了
#===============================================

#==========================================
# &Subroutines
#==========================================
sub ModifyeHtunerInput
{
	$App->print("\n\nModify eHtuner input:\n");

	my $InFile     = $Args->GetGetArg(0);
	$InFile = 'eHparam.ini' if(!$InFile);
	my $OutFile    = "$InFile.modified";
	my $BindPath   = $Args->GetGetArg(1);
	my $OS = $^O;
	if(!$BindPath) {
		if($OS eq 'MSWin32') {
			$BindPath = 'bind.exe';
		}
		elsif($OS eq 'linux') {
			$BindPath = 'bind';
		}
	}
	my $FitDOSPath = $Args->GetGetArg(2);
	if(!$FitDOSPath) {
		if($OS eq 'MSWin32') {
			$FitDOSPath = 'fit_dos.exe';
		}
		elsif($OS eq 'linux') {
			$FitDOSPath = 'fit_dos';
		}
	}

	$App->print("Input (eHtuner created): $InFile\n");
	$App->print("Output file            : $OutFile\n");
	$App->print("bind path              : $BindPath\n");
	$App->print("fit_dos path           : $FitDOSPath\n");

	my $in  = JFile->new($InFile, 'rb')  or die "$!: Can not read [$InFile].\n";
	my $out = JFile->new($OutFile, 'wb') or die "$!: Can not write to [$OutFile].\n";
	while(1) {
		my $line = $in->ReadLine();
		last if(!defined $line);
		
		if($line =~ /^BIND_COMMAND/i) {
			$out->print($line);
			$in->ReadLine();
			$out->print("$BindPath\n");
		}
		elsif($line =~ /^FITDOS_COMMAND/i) {
			$out->print($line);
			$in->ReadLine();
			$out->print("$FitDOSPath\n");
		}
		else {
			$out->print($line);
		}
	}
	$out->Close();
	$in->Close();

	$App->print("\n\nModify eHtuner input: finished\n");
}

sub MakeDOSCSV
{
	$App->print("\n\nMake DOS CSV from output of YAeHMOP:\n");

	my $InFile = $Args->GetGetArg(0);
	$InFile = 'bind.in.out' if(!$InFile);
	my $CSVFile = $Args->GetGetArg(1);
	$CSVFile = 'DOS.csv' if(!$CSVFile);

	$App->print("Input (bind output): $InFile\n");
	$App->print("DOS CSV file       : $CSVFile\n");

	my $in  = JFile->new($InFile, 'r') or die "$!: Can not read [$InFile]\,";
# ### TOTAL DENSITY OF STATES 
# 16 states are present.
# #BEGIN CURVE
#  0.096000000 -27.502025604
# #END CURVE
#   0.040943685 -27.502025604
# #END CURVE
	$in->SkipTo("### TOTAL DENSITY OF STATES");
	$in->ReadLine();
	$in->ReadLine();
	my $c = 0;
	my @DOS;
	while(1) {
		my $line = $in->ReadLine();
		last if(!defined $line or $line =~ /^\s*#\s*END/i);

		my @a = Utils::Split("\\s+", $line);
		last if(@a < 2);

		$DOS[0][$c++] = [$a[1], $a[0]];
#print "c1=$c\n";
	}
# ### PROJECTED DENSITY OF STATES 
# 2.000000 states are present.
# ORBITAL 1 1.000000,
# #BEGIN CURVE
	my $id  = 1;
	my $End = 0;
	while(1) {
		$in->SkipTo("### PROJECTED DENSITY OF STATES");
		$in->ReadLine();
		$in->ReadLine();
		$in->ReadLine();
		$c = 0;
		while(1) {
			my $line = $in->ReadLine();
			if(!defined $line or $line =~ /# END OF DOS/i) {
				$End = 1;
				last;
			}
			last if($line =~ /^\s*#\s*END/i);

			my @a = Utils::Split("\\s+", $line);
			last if(@a < 2);

			$DOS[$id][$c++] = [$a[1], $a[0]];
#print "$id c2=$c\n";
		}
		$id++;
		last if($End);
	}
	$in->Close();

	my $nData = @{$DOS[0]};
print "nData=$nData\n";	
	my $out = JFile->new($CSVFile, 'w') or die "$!: Can not write to [$CSVFile]\,";
	$out->print("E(eV),Total DOS");
	for(my $id = 1 ; $id < @DOS ; $id++) {
		$out->print(",E(eV),DOS $id");
	}
	$out->print("\n");
	for(my $i = 0 ; $i < $nData ; $i++) {
		for(my $id = 0 ; $id < @DOS ; $id++) {
			my $p = $DOS[$id];
			$out->print("$p->[$i][0],$p->[$i][1],");
#print "$id: $i: ($p->[$i][0], $p->[$i][1])\n";
		}
		$out->print("\n");
	}
	$out->Close();

	$App->print("\n\nMake DOS CSV from output of YAeHMOP: finiahed\n");
}

sub MakeInput
{
	$App->print("\n\nMake YAeHMOP input File:\n");

	my $aKProduct = $Args->GetGetArg('aKProduct', 1, 1);
#	my $aKProduct = $Args->GetGetArg('aKProduct', 0);
#	$aKProduct = 2.0 if(!$aKProduct);
	my $eHParamdb = $Args->GetGetArg('db',        1, 1);
#	my $eHParamdb = $Args->GetGetArg('db',        0);
	$eHParamdb = 'eht_parms.dat' if(!$eHParamdb);
	my $KPOINTS   = $Args->GetGetArg('KPOINTS',   0);
	if(!$KPOINTS) {
		my @files = (glob("*.klist"), glob("*.KPOINTS"));
		$KPOINTS = $files[0] if(@files);
	}
	my $IBZKPT    = $Args->GetGetArg('IBZKPT',    1, 1);
#	my $IBZKPT    = $Args->GetGetArg('IBZKPT',    0);
	$IBZKPT = 'IBZKPT' if(!$IBZKPT);
	my $CIFFile   = $Args->GetGetArg(0);
	if($CIFFile eq 'auto') {
		my @files = (glob('*.cif'), glob('*.CIF'));
		$CIFFile = $files[0] if(@files);
	}
	my $InpFile   = $Args->GetGetArg('Output');
	$InpFile = 'bind.in' if(!$InpFile);

	unless($CIFFile) {
		$App->print("File name should be specified.\n");
		$App->Usage();
		return 0;
	}

	if(!-f $CIFFile) {
		$CIFFile = '*.cif' if($CIFFile eq '');
		my @f = glob($CIFFile);
		$CIFFile = $f[0];
	}
#ファイル名を（ベース名, ディレクトリ名, 拡張子）に分解
	my ($drive, $dir, $filename, $ext, $lastdir, $filebody) = Deps::SplitFilePath($CIFFile);
	my $SampleName = $filebody;
	$InpFile = "$SampleName.in" if(!defined $InpFile);
	$App->print("CIF File     : $CIFFile\n");
	$App->print("Output       : $InpFile\n");
	$App->print("eHParam db   : $eHParamdb\n");
	$App->print("Band k points: $KPOINTS\n");
	$App->print("SCF  kpoints : $IBZKPT\n");
	$App->print("aKProduct    : $aKProduct\n");

	$App->print("  Read [$CIFFile].\n");

	my $CIF = new CIF();
	unless($CIF->Read($CIFFile)) {
		$App->print("Error: Can not read $CIFFile.\n\n");
		return 0;
	}

#CIFクラスの内容から、Crystalクラスを作成
	my $Crystal = $CIF->GetCCrystal();
	$Crystal->ExpandCoordinates();
#($this, $LatticeParameterOnly, $CheckConsistency, $IsPrint) = @_;
	my $PrimCrystal = $Crystal->GetPrimitiveCrystal();
	my $T           = $Crystal->GetMatrixConventionalToPrimitiveCell(0);

	my ($a,$b,$c,$alpha,$beta,$gamma) = $Crystal->LatticeParameters();
#	my @AtomTypeList = $Crystal->GetCAtomTypeList();
#	my @SiteList     = $Crystal->GetCAsymmetricAtomSiteList();
#	my @SiteList     = $Crystal->GetCExpandedAtomSiteListByOutputMode();
#	my @SiteList     = $Crystal->GetCExpandedAtomSiteList();
	my @AtomTypeList = $PrimCrystal->GetCAtomTypeList();
	my @SiteList     = $PrimCrystal->GetCExpandedAtomSiteList();
	my $nSite = @SiteList;

	my $band = new BandStructure;
	my @klist;
	if($KPOINTS) {
		$App->print("  Read [$KPOINTS].\n");
		@klist = $band->ReadKListFile($KPOINTS);
	}
	if(@klist == 0) {
		$App->print("  KPOINTS [$KPOINTS] is not valid. Use default:\n");
		push(@klist, ['GAMMA', 0.0, 0.0, 0.0]);
		push(@klist, ['X',     0.5, 0.0, 0.0]);
		push(@klist, ['R',     0.5, 0.5, 0.5]);
		push(@klist, ['M',     0.5, 0.5, 0.0]);
	}
	my $nK = @klist;
for(my $i = 0 ; $i < @klist ; $i++) {
print "    $i: $klist[$i]{name}: $klist[$i]{kx}, $klist[$i]{ky}, $klist[$i]{kz}\n";
}

	my @SCFklist;
	if($IBZKPT) {
		$App->print("  Read [$IBZKPT].\n");
		my $in = JFile->new($IBZKPT, 'r');
		if($in) {
			$in->SkipTo('Reciprocal lattice');
			my $c = 0;
			while(1) {
				my $line = $in->ReadLine();
				my @aa = Utils::Split("\\s+", $line);
				last if(@aa < 4);
				$SCFklist[$c] = [$aa[0], $aa[1], $aa[2], $aa[3]];
				$c++;
			}
			$in->Close();
		}
	}
	my $nSCFK = @SCFklist;
for(my $i = 0 ; $i < @SCFklist ; $i++) {
print "    $i: $SCFklist[$i][0], $SCFklist[$i][1], $SCFklist[$i][2], $SCFklist[$i][3]\n";
}

	my ($nx, $ny, $nz, $NKREDX, $NKREDY, $NKREDZ) 
		= $band->AnalyzeaKProduct($Crystal, $aKProduct);

	$App->print("  Save to [$InpFile].\n");
	my $out = JFile->new($InpFile, 'wb') or die "$!: Can not write to [$InpFile].\n";
	$out->print(";------------------------------------------------------------------------\n");
	$out->print("\n");
	$out->print("$SampleName\n");
	$out->print("\n");
	$out->print("; Define the geometry in crystallographic coordinates\n");
	$out->print("Geometry Crystallographic\n");
	$out->print("; Number of atoms\n");
	$out->printf("%d\n", $nSite+4);
	$out->print("; '*' indicates that custom parameters are being used.\n");
	$out->print("; '&' are dummy coordinates to define the unit cell.\n");
#$out->print("1 * 0.00000  0.00000  0.0000\n");
	my %aname;
	my @aname;
	my %eHParam;
	for(my $i = 0 ; $i < $nSite ; $i++) {
		my $atom        = $SiteList[$i];
		my $name        = $atom->AtomNameOnly();
		my ($x, $y, $z) = $atom->Position();
		$eHParam{uc $name} = &ReadeHParamDB($eHParamdb, $name) if(!defined $eHParam{uc $name});
		my $idx = $i+1;
		if($aname{$name}) {
			$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, $name, $x, $y, $z);
		}
		else {
			$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '*', $x, $y, $z);
			push(@aname, $name);
		}
		$aname{$name}++;
	}
	my $idx = $nSite + 1;
	$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', 0.0, 0.0, 0.0);
	if($T) {
		$idx++;
		$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', $T->[0][0], $T->[0][1], $T->[0][2]);
		$idx++;
		$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', $T->[1][0], $T->[1][1], $T->[1][2]);
		$idx++;
		$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', $T->[2][0], $T->[2][1], $T->[2][2]);
	}
	else {
		$idx++;
		$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', 1.0, 0.0, 0.0);
		$idx++;
		$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', 0.0, 1.0, 0.0);
		$idx++;
		$out->printf("%2d %2s %7.5f  %7.5f  %7.5f\n", $idx, '&', 0.0, 0.0, 1.0);
	}

	$out->print("\n");
	$out->print("; Define custom eH parameters \n");
	$out->print("; Symbol AtomNum nElectrons [s] [p] [d] [f]\n");
	$out->print("; [s], [p]: n zeta IP\n");
	$out->print("; [d], [f]: n zeta1 IP c1 zeta2 c2\n");
	$out->print("parameters\n");
#for my $key (keys %eHParam) {
#print "key=$key\n";
#}
	my $nTotalVEL = 0.0;
	for(my $i = 0 ; $i < $nSite ; $i++) {
		my $atom        = $SiteList[$i];
		my $name        = $atom->AtomNameOnly();
		my $p    = $eHParam{uc $name};
		$nTotalVEL += $p->{s}{NVElectron};
#print "i=$i: $atom: ne=$p->{s}{NVElectron}; tot=$nTotalVEL\n";
	}
	for(my $i = 0 ; $i < @aname ; $i++) {
		my $name = $aname[$i];
		my $p    = $eHParam{uc $name};
		my $nVEL = $p->{s}{NVElectron};
#		$nVEL += $p->{p}{NVElectron} if($p->{p});
#		$nVEL += $p->{d}{NVElectron} if($p->{d});
#		$nVEL += $p->{f}{NVElectron} if($p->{f});
#print "AT=$name   p=$p  ps=$p->{s} nVEL=$nVEL\n";
		$out->printf("%2s %3d %3d ", $name, $p->{s}{AtomNum}, $nVEL);
#		$out->printf("%2s %3d %6.4f ", $name, $p->{s}{AtomNum}, $nVEL);
		if($p->{p}) {
			$out->printf(" %2d %8.4f %8.4f ", $p->{s}{nQuantum}, $p->{s}{exp1}, $p->{s}{IP});
		}
		if($p->{p}) {
			$out->printf(" %2d %8.4f %8.4f ", $p->{p}{nQuantum}, $p->{p}{exp1}, $p->{p}{IP});
		}
		if($p->{d}) {
			$out->printf(" %2d %8.4f %8.4f %8.4f %8.4f %8.4f ", 
				$p->{d}{nQuantum}, $p->{d}{exp1}, $p->{d}{IP}, $p->{d}{coeff1}, 
				$p->{d}{exp2}, $p->{d}{coeff2});
		}
		if($p->{f}) {
			$out->printf(" %2d %8.4f %8.4f %8.4f %8.4f %8.4f ", 
				$p->{f}{nQuantum}, $p->{f}{exp1}, $p->{f}{IP}, $p->{f}{coeff1}, 
				$p->{f}{exp2}, $p->{f}{coeff2});
		}
		$out->print("\n");
	}	
	$out->print("\n");
	$out->print("; Define the lattice\n");
	$out->print("lattice\n");
	$out->print("; Dimensionality\n");
	$out->print("3\n");
	$out->print("; Number of overlaps in each direction (~1000 atoms size is recommended)\n");
	$out->print("$nx $ny $nz\n");
	$out->print("; definition of a\n");
	$out->printf("%d %d\n", $nSite+1, $nSite+2);
	$out->print("; definition of b\n");
	$out->printf("%d %d\n", $nSite+1, $nSite+3);
	$out->print("; definition of c\n");
	$out->printf("%d %d\n", $nSite+1, $nSite+4);
	$out->print("\n");
	$out->print("Crystal Spec\n");
	$out->print("; Lengths of the lattice vectors a, b, and c\n");
	$out->print("$a $b $c\n");
	$out->print("; Lattice angles alpha, beta, and gamma\n");
	$out->print("$alpha $beta $gamma\n");
	$out->print("\n");
	$out->print("Band\n");
	$out->print("40\n");
	$out->print("$nK\n");
	$out->print("; The k points are specified by the reciprocal lattice of the primitive cell.\n");
	$out->print("; The conversions between the two cells can be found in Ashcroft and Mermin\n");
	for(my $i = 0 ; $i < $nK ; $i++) {
		$out->printf("%6s  %8.4f %8.4f %8.4f\n", 
			$klist[$i]{name}, $klist[$i]{kx}, $klist[$i]{ky}, $klist[$i]{kz});
	}
	$out->print("\n");
	$out->print("print\n");
	$out->print("distance matrix\n");
	$out->print("overlap matrix\n");
	$out->print("hamiltonian\n");
	$out->print("wave func transp\n");
	$out->print("charge mat transp\n");
	$out->print("energy levels\n");
	$out->print("end_print\n");
	$out->print("\n");
	if($nSCFK > 0) {
		$out->print("average properties\n");
	}
	else {
		$out->print(";average properties\n");
	}
	$out->print("\n");
	$out->print("electrons\n");
	$out->print("$nTotalVEL\n");
	$out->print(";charge\n");
	$out->print(";0\n");
	$out->print("\n");

	$out->print("; write output that cooperate can use\n");
	$out->print(";dump distance matrix\n");
	$out->print("\n");
	my $nAtomType = @aname;
	$out->print("projected dos\n");
	$out->print("$nAtomType\n");
	my $c = 1;
	for(my $i = 0 ; $i < @aname ; $i++) {
		$out->print(";$aname[$i] ");
		my $p = $eHParam{uc $aname[$i]};
		$out->printf("%d:%s ", $c++, 's') if($p->{s});
		$out->printf("%d:%s %d:%s %d:%s ", $c++, 'px', $c++, 'py', $c++, 'pz') if($p->{p});
		$out->printf("%d:%s %d:%s %d:%s %d:%s %d:%s %d:%s %d:%s "
				, $c++, 'dx2-y2', $c++, 'dz2', $c++, 'dxy', $c++, 'dxz', $c++, 'dyz') if($p->{d});
		$out->printf("%d:%s %d:%s %d:%s %d:%s %d:%s %d:%s %d:%s %d:%s %d:%s "
				, $c++, 'f1', $c++, 'f2', $c++, 'f3', $c++, 'f4', $c++, 'f5', $c++, 'f6', $c++, 'f7') if($p->{f});
		$out->print("\n");
	}
	$c = 1;
	for(my $i = 0 ; $i < @aname ; $i++) {
		$out->print("Orbital ");
		my $p = $eHParam{uc $aname[$i]};
		for(my $j = 0 ; $j < $nSite ; $j++) {
			my $atom        = $SiteList[$j];
			my $name        = $atom->AtomNameOnly();
			if(uc $name eq uc $aname[$i]) {
				$out->printf("%d %3.1f  ", $c++, 1.0) if($p->{s});
				$out->printf("%d %3.1f %d %3.1f %d %3.1f  ", $c++, 1.0, $c++, 1.0, $c++, 1.0) if($p->{p});
				$out->printf("%d %3.1f %d %3.1f %d %3.1f %d %3.1f %d %3.1f  "
					, $c++, 1.0, $c++, 1.0, $c++, 1.0, $c++, 1.0, $c++, 1.0) if($p->{d});
				$out->printf("%d %3.1f %d %3.1f %d %3.1f %d %3.1f %d %3.1f %d %3.1f %d %3.1f  "
					, $c++, 1.0, $c++, 1.0, $c++, 1.0, $c++, 1.0, $c++, 1.0, $c++, 1.0, $c++, 1.0) if($p->{f});
			}
		}
		$out->print("\n");
	}
	for(my $i = 0 ; $i < @aname ; $i++) {
		$out->print(";atom ");
		for(my $j = 0 ; $j < $nSite ; $j++) {
			my $atom        = $SiteList[$j];
			my $name        = $atom->AtomNameOnly();
			if(uc $name eq uc $aname[$i]) {
				$out->printf(" %2d 1.0 ", $j+1);
			}
		}
		$out->print("\n");
	}
	$out->print("\n");
	$out->print("; These k points are for average properties calculations\n");
	$out->print("; kx ky kz weight\n");
	if($nSCFK > 0) {
		$out->print("k points\n");
		$out->print("$nSCFK\n");
		for(my $i = 0 ; $i < $nSCFK ; $i++) {
			$out->printf("%7.4f %7.4f %7.4f  %2d\n",
				$SCFklist[$i][0], $SCFklist[$i][1], $SCFklist[$i][2], $SCFklist[$i][3]);
		}
	}
	else {
		$out->print(";k points\n");
		$out->print(";32\n");
		$out->print(";-.33 -.33 .33 2\n");
	}
	$out->print("\n");

	$out->Close();

	$App->print("\n***Make YAeHMOP input file: Finished.\n\n");

	return;
}

#;EHParameters
#;numberOfBasisFunctions: 300
#;Atomlab AtNo Nvalen Nzeta Nquant Ang  IP       exp1     exp2   coeff1   coeff2  
#CE      58     4       1     6     s   -4.970   1.7990   0.0000   1.0000   0.0000 
#CE      58     4       1     6     p   -4.970   1.7990   0.0000   1.0000   0.0000 
#CE      58     4       1     5     d   -8.430   2.7470   0.0000   1.0000   0.0000
#CE      58     4       1     4     f  -11.280   3.9070   0.0000   1.0000   0.0000 

#Symbol AtomNum NumValEl ns zetas IPs np zetap IPp nd zeta1d IPd c1 zeta2d c2
#   f is the same as d
#Ti 22 4 4 1.075 -8.970 4 1.075 -5.400 3 4.55 -10.81 .4206 1.400 .7839
sub ReadeHParamDB
{
	my ($db, $atomname) = @_;
	
	print "  Read [$db].\n";

	my $in = JFile->new($db, 'r') or die "$!: Can not read [$db].\n";

	my $an = uc $atomname;
print "atom: $an\n";
	my %Hash;
	while(1) {
		my $line = $in->ReadLine();
		last if(!defined $line);
		next if($line =~ /^\s*;/);

		my @aa = Utils::Split("\\s+", $line);
		next if($an ne uc $aa[0]);
print "   Hit: $line";
		$Hash{$aa[5]} = {
			Name       => $aa[0],
			AtomNum    => $aa[1],
			NVElectron => $aa[2],
			nZeta      => $aa[3],
			nQuantum   => $aa[4],
			Angular    => $aa[5],
			IP         => $aa[6],
			exp1       => $aa[7],
			exp2       => $aa[8],
			coeff1     => $aa[9],
			coeff2     => $aa[10],
			};
	}
	
	$in->Close();
	
	return \%Hash;
}