#!/usr/bin/perl

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

use strict;

use Utils;
use MyApplication;
use Crystal::WIEN2k;
use Crystal::XCrySDen;

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

#===============================================
# スクリプト大域変数
#===============================================
# for klist & bxsf
my $div = 360;
my ($nx, $ny, $nz)  = (20, 20, 40);
my $Target = "RE"; # RE|IM|ABS|ARG|PSI
my $iKPoint = 1;
my $iBand   = 0;

#==========================================
# コマンドラインオプション読み込み
#==========================================
$App->AddArgument("--Action",    "--Action=[all]", '');
$App->AddArgument("--DebugMode", "--DebugMode: Set DebugMode", '');
exit 1 if($App->ReadArgs(0) != 1);
my %Args = $App->Args()->GetArgHash();

#==========================================
# メイン関数スタート
#==========================================
my $Debug = $Args{DebugMode};
$App->SetDebug($Debug);
my $Action  = $Args{Action};
my $DataDir = $App->Args()->GetGetArg(0);
$DataDir = "." if($DataDir eq '');
$DataDir = Utils::ReduceDirectory($DataDir);
my $StructFile = &GetStructPath($DataDir);
if(!-e $StructFile) {
	$App->print("Error: StructFile [$StructFile] could not obtained.\n");
	exit 0;
}
my $Title    = &GetTitle($StructFile);
my $In7Path  = Deps::ReplaceExtension($StructFile, ".in7");
my $RhoPath  = Deps::ReplaceExtension($StructFile, ".psink");
my $VASPPath = Deps::ReplaceExtension($StructFile, ".${Target}.vasp");
my $PrefPath = Deps::MakePath($DataDir, ".prefs-scf");

my $WIEN2k = new WIEN2k;
my $Crystal = $WIEN2k->ReadStructFile($StructFile, 0);
unless($Crystal) {
	$App->print("Error: Can not read [$StructFile].\n");
	return 0;
}

my $config = new IniFile();
$config->ReadAll($PrefPath);

$App->print("Struct: $StructFile\n");
$App->print("Lattice System: ", $Crystal->LatticeSystem(), "\n");
$App->print("Space Group   : ", $Crystal->SPGName(), "\n");
if($config->{spin} eq '') {
	$App->print("Spin polarized: no\n");
}
else {
	$App->print("Spin polarized: yes\n");
}

my $ret = 0;
if($Action eq 'all') {
	&Execute($Crystal, $config, $Title, $Target, $VASPPath, $In7Path, $RhoPath);
}
else {
	$App->print("Error: Invald Action: [$Action]\n");
	$App->Usage();
}

exit $ret;

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

#==========================================
# &Subroutines
#==========================================
sub GetTitle
{
	my ($path) = @_;

	my ($drive, $directory, $filename, $ext, $lastdir, $filebody) = Deps::SplitFilePath($path);
	my $Title   = $filebody;
	$Title =~ s/\s/_/g;
	return $Title;
}

sub GetStructPath
{
	my ($DataDir) = @_;

	my ($drive, $directory, $filename, $ext, $lastdir, $filebody) = Deps::SplitFilePath($DataDir);
	my $StructFile = Deps::MakePath($DataDir, "$filename.struct", 0);
	if(!-e $StructFile) {
		my $fmask = Deps::MakePath($DataDir, "*.struct");
		my @files = glob($fmask);
		if($files[0] =~ /setrmt\.struct$/) {
			$StructFile = $files[1];
		}
		else {
			$StructFile = $files[0];
		}
		unless($StructFile) {
			print("Struct file is not found from Data Directory [$DataDir].\n");
			return undef;
		}
	}
	return $StructFile;
}

sub Execute
{
	my ($Crystal, $config, $Title, $Target, $VASPPath, $In7Path, $RhoPath) = @_;

	$App->print("\nMake Wavefunction .vasp file\n");
	$App->print("VASP Path: $VASPPath\n");

	unlink($VASPPath);
	my $nx1 = $nx+1;
	my $ny1 = $ny+1;
	my $nz1 = $nz+1;
	my @Rho;
	for(my $iz = 0 ; $iz <= $nz ; $iz++) {
$App->print("iz=$iz / $nz\n");
		unlink($RhoPath);
		open(IN5, ">$In7Path") or die "Error: $!: Can not write to [$In7Path].\n";

		my $iz1 = int($iz * $div / $nz + 0.01);
		print IN5 <<EOT;
2D NON-ORTHOGONAL
   0    0  $iz1  $div   # origin of plot: x, y,z, denominator
 $div    0  $iz1  $div   # x-end of plot
   0  $div  $iz1  $div   # y-end of plot
 $nx1  $ny1  1  1       # nx, ny, xincrement, yincrement
EOT
		printf(IN5 "%-s                    # DEP(HASING)|NO(POST-PROCESSING)\n", "NO");
		printf(IN5 "%-3s %-3s %-5s           # switch ANG|ATU|AU  LARGE|SMALL\n", $Target, "ANG", "LARGE");
		print(IN5 "$iKPoint $iBand         # k-point, band index\n");
		close(IN5);
		
		my $cmd = ($config->{spin} eq '')? "x lapw7" : "x lapw7 -up";
		if(Utils::Execute($cmd)) {
			$App->print("Error: Could not execute [$cmd].\n");
			exit 0;
		}
		
		my $in = new JFile();
		if(!$in->Open($RhoPath, "r")) {
			$App->print("Error: Could not read [$RhoPath].\n");
			exit 0;
		}
		my $line = $in->ReadLine();
		while(1) {
			$line = $in->ReadLine();
			last if(!defined $line);

			push(@Rho, Utils::Split("\\s+", $line))
		}
		$in->Close();
	}
	unlink($RhoPath);

	my $out = new JFile();
	if(!$out->Open($VASPPath, "w")) {
		$App->print("Error: Could not write to [$VASPPath].\n");
		exit 0;
	}

	my ($a11, $a12, $a13, $a21, $a22, $a23, $a31, $a32, $a33) = $Crystal->LatticeVectors();
	$out->print("$Title\n");
	$out->printf("%6.4f\n", $a11);
	$out->printf("   %10.6f  %10.6f  %10.6f\n", 1.0,       $a12/$a11, $a13/$a11);
	$out->printf("   %10.6f  %10.6f  %10.6f\n", $a21/$a11, $a22/$a11, $a23/$a11);
	$out->printf("   %10.6f  %10.6f  %10.6f\n", $a31/$a11, $a32/$a11, $a33/$a11);
	$out->print("    1\n");
	$out->print("Direct\n");
	$out->print("  0.000000  0.000000  0.000000\n");
	$out->print("\n");
	$out->printf(" %d %d %d\n", $nx1, $ny1, $nz1);
	my $count = 0;
	for(my $ix = 0 ; $ix <= $nx ; $ix++) {
		for(my $iy = 0 ; $iy <= $ny ; $iy++) {
			for(my $iz = 0 ; $iz <= $nz ; $iz++) {
				$out->printf("%8.5f    ", $Rho[$count]);
				$count++;
				if($count % 10 == 0) {
					$out->print("\n");
				}
			}
		}
	}
	$out->Close();

}