#!/usr/bin/perl

use lib 'd:/Programs/Perl/lib';

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

use strict;
#use warnings;

use Math::Matrix;
use Math::MatrixReal;

use Deps;
use Utils;
use JFile;

use MyApplication;
use CSV;
use Sci::Science;
use Sci::Algorism;

#===============================================
# デバッグ関係変数
#===============================================
#$PrintLevelが大きいほど、情報が詳しくなる
my $PrintLevel = 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();

my $nDivide = 1;
my $ExtraporateLow  = "Drude";
#my $ExtraporateLow  = "None";
my $ExtraporateHigh = "Constant";
#my $ExtraporateHigh = "None";
#my $Integrator = "Rectangular";
my $Integrator = "Simpson";

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

#==========================================
# メイン関数スタート
#==========================================
#Utils::InitHTML("Research", $WebCharSet, "_self");

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

my $ret = 0;

&Execute();

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

#Utils::EndHTML();

exit $ret;

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

#==========================================
# &Subroutines
#==========================================
sub Execute
{
	my $h  = Sci::h();
	my $e  = Sci::e();
	my $pi = Sci::pi();

	my $InputFile  = $Args->GetGetArg(0);
	$InputFile = "Au.asp.csv" if(!defined $InputFile);
	my $OutputFile = $Args->GetGetArg(1);
	if(!defined $OutputFile or $InputFile eq $OutputFile) {
		$OutputFile = Deps::ReplaceExtension($InputFile, "KK.csv");
	}

	$App->print("  Input file : $InputFile\n");
	$App->print("  Output file: $OutputFile\n");
	$App->print("  ExtraporateLow : $ExtraporateLow\n");
	$App->print("  ExtraporateHigh: $ExtraporateHigh\n");
	$App->print("  Integrator     : $Integrator\n");

	my $csv = new CSV;
	if($csv->Read($InputFile)) {
		$App->print("  Read [$InputFile]\n");
	}
	else {
		$App->print("  Error!!: Can not read [$InputFile]\n");
		return -1;
	}

	my $pLabelArray = $csv->LabelArray();
	my $pDataArray  = $csv->DataArray();
	my $nDataRow = @$pDataArray;
	my $pX       = $pDataArray->[0];
	my $nData = @$pX;
	$App->print("  nDataRow: $nDataRow\n");
	$App->print("  nData   : $nData\n");
	my $x0 = $pX->[0];
	my $x1 = $pX->[@$pX - 1];
	my $dx = $pX->[1] - $pX->[0];
	$App->print("  x: $x0 - $x1     dx: $dx\n");
#	for(my $i = 0 ; $i < $nDataRow ; $i++) {
#		$App->print("  Row $i: ", $pLabelArray->[$i], "\n");
#	}

	my $Eidx = -1;
	for(my $i = 0 ; $i < $nDataRow ; $i++) {
		if($pLabelArray->[$i] eq 'E' or $pLabelArray->[$i] =~ /^E[ \(]*eV[ \)]*/ ) {
			$Eidx = $i;
			last;
		}
	}
	$App->print("  E data is at Row $Eidx\n");
	my $Ridx = -1;
	for(my $i = 0 ; $i < $nDataRow ; $i++) {
		if($pLabelArray->[$i] eq 'R') {
			$Ridx = $i;
			last;
		}
	}
	$App->print("  R data is at Row $Ridx\n");

	my $nidx = -1;
	for(my $i = 0 ; $i < $nDataRow ; $i++) {
		if($pLabelArray->[$i] eq 'n') {
			$nidx = $i;
			last;
		}
	}
	$App->print("  n data is at Row $nidx\n") if($nidx >= 0);
	my $kidx = -1;
	for(my $i = 0 ; $i < $nDataRow ; $i++) {
		if($pLabelArray->[$i] eq 'k') {
			$kidx = $i;
			last;
		}
	}
	$App->print("  k data is at Row $kidx\n") if($kidx >= 0);

	my $pf = $pDataArray->[$Eidx];
	my $pR = $pDataArray->[$Ridx];
	for(my $i = 0 ; $i < $nData ; $i++) {
		$pf->[$i] = $pf->[$i] * $e / $h;
	}

	my $Minf = $pf->[0];
	my $Maxf = $pf->[$nData-1];
	my $df = ($Maxf - $Minf) / ($nData - 1);
	$App->printf("  f range: %8.4e - %8.4e Hz (%8.4e Hz)\n", $Minf, $Maxf, $df);

	my $nDataOrg = $nData;
	$nData *= $nDivide;
print " nData: $nDataOrg => $nData\n";
	my $dfOrg = $df;
	$df /= $nDivide;
print " delta f: $dfOrg => $df\n";
	my @f;
	my @R;
	my @lnR2;
	my @nOrgInt;
	my @kOrgInt;
	for(my $i = 0; $i < $nData ; $i++) {
		my $x = $Minf + $i * $df;
		my $v = Algorism::Interporate($nDataOrg, $Minf, $dfOrg, $pR, $x);
#print "$i: ($x, $v)\n";
		$f[$i] = $x;
		$v = 1e-10 if($v <= 1e-10);
		$R[$i] = $v;
#my $E = $h * $x / $e;
#print "$i: $E   $x  $R[$i]\n";
		$lnR2[$i] = 0.5 * log($v);
		if($nidx >= 0) {
			$nOrgInt[$i] = Algorism::Interporate(
					$nDataOrg, $Minf, $dfOrg, $pDataArray->[$nidx], $x);
		}
		if($kidx >= 0) {
			$kOrgInt[$i] = Algorism::Interporate(
					$nDataOrg, $Minf, $dfOrg, $pDataArray->[$kidx], $x);
		}
	}

	my @Q;
	my @n;
	my @k;
	my $Coeff = 2.0/$pi * $df;
	for(my $i = 0 ; $i < $nData ; $i++) {
		my $fg = $f[$i];
		my $idx = int(($fg - $Minf) / $df + 0.001);
		my $j0 = 0;
		$j0 = 1  if($idx % 2 == 0);

		my $v = 0.0;
		if($Integrator eq "Rectangular") {
			for(my $j = $j0 ; $j < $nData ; $j += 2) {
				my $fj = $f[$j];
#				$v += 2.0 * $lnR2[$j] / ($fj*$fj - $fg*$fg);
				$v += log($R[$j]/$R[$i]) / ($fj*$fj - $fg*$fg);
#print "$j: $v\n" if($i == 100);
			}
#return 0 if($i == 100);
		}
		elsif($Integrator eq "Simpson") {
			for(my $j = $j0 ; $j <= $nData-3 ; $j += 2) {
				next if($i-3 <= $j and $j <= $i);
				my $fj = $f[$j];
				my $fj1 = $f[$j+1];
				my $fj2 = $f[$j+2];
				$v += 1.0/3.0 * (
					0.5 * log($R[$j]/$R[$i])   / ($fj*$fj   - $fg*$fg)
				      + 2.0 * log($R[$j+1]/$R[$i]) / ($fj1*$fj1 - $fg*$fg)
				      + 0.5 * log($R[$j+2]/$R[$i]) / ($fj2*$fj2 - $fg*$fg)
					);
#				$v += 1.0/3.0 * (
#					          $lnR2[$j]   / ($fj*$fj   - $fg*$fg)
#				      + 4.0 * $lnR2[$j+1] / ($fj1*$fj1 - $fg*$fg)
#				      +       $lnR2[$j+2] / ($fj2*$fj2 - $fg*$fg)
#					);
			}
		}
		else {
			$App->print("Error!!: Invalid Integrator [$Integrator]\n");
			return -5;
		}
		$v *= $Coeff * $fg;

		my $jmax = $nData - 1;
		$jmax-- if(($j0+$jmax) % 2 != 0);
		$j0 = 0 if($j0 == 1);
		my $fmin    = $f[$j0] - $df;
		$fmin -= $df if(abs($fmin-$fg) < 1.0);
		my $lnR2min = $lnR2[$j0];
		my $fmax    = $f[$jmax] + $df;
		$fmax += $df if(abs($fmax-$fg) < 1.0);
		my $lnR2max = $lnR2[$nData-1];
#		my $fmin    = $f[0] - $df;
#		my $fmax    = $f[$nData-1] + $df;
#		my $lnR2min = $lnR2[0];
#		my $lnR2max = $lnR2[$nData-1];

		if($ExtraporateLow  eq "Drude") {
			$v += $lnR2min / $pi * log( ($fg-$fmin)/($fmin+$fg) )
		}
		elsif($ExtraporateLow  eq "None") {
		}
		else {
			$App->print("Error!!: Invalid Extraporation to low f [$ExtraporateLow]\n");
			return -3;
		}
		if($ExtraporateHigh eq "Constant") {
			$v += $lnR2max / $pi * log( ($fmax+$fg)/($fmax-$fg) );
		}
		elsif($ExtraporateHigh  eq "None") {
		}
		else {
			$App->print("Error!!: Invalid Extraporation to high f [$ExtraporateHigh]\n");
			return -3;
		}

		$Q[$i] = -$v;
		my $SqrR = sqrt($R[$i]);
		my $div = 1.0 + $R[$i] - 2.0 * $SqrR * cos($Q[$i]);
		$n[$i] = (1.0 - $R[$i]) / $div;
		$k[$i] = 2.0 * $SqrR * sin($Q[$i]) / $div;
	}

	my @LabelArray = ("f(s-1)", "R", "0.5lnR", "Q", "n", "k");
	if($nidx >= 0) {
		push(@LabelArray, "n(org)");
	}
	if($kidx >= 0) {
		push(@LabelArray, "k(org)");
	}
	my @DataArray  = (\@f, \@R, \@lnR2, \@Q, \@n, \@k);
	if($nidx >= 0) {
		push(@DataArray, \@nOrgInt);
	}
	if($kidx >= 0) {
		push(@DataArray, \@kOrgInt);
	}

	$csv->SetLabelArray(\@LabelArray);
	$csv->SetDataArray(\@DataArray);
	if($csv->Save($OutputFile)) {
		$App->print("  Save to [$OutputFile]\n");
	}
	else {
		$App->print("  Error!!: Can not write to [$OutputFile]\n");
		return -2;
	}

	$App->print("\n\n<b>Convolute CSV File: Finished</b>\n");
}