#!/usr/bin/perl

use lib 'd:/Programs/Perl/lib';
#use lib 'd:/MyWebs/cgi-bin/lib';

use lib '/home/tkamiya/bin/lib';
use lib '/home/tkamiya/bin';

use strict;
#use warnings;
use Math::MatrixReal;
use Math::Complex;

use Sci qw($pi $a0 $torad $todeg);# acos asin);
use Sci::Algorism;
use Sci::MyMatrixReal;

use MyApplication;
use Crystal::CIF;
use Crystal::Crystal;

#No P1 RhombHex HexRhomb HexOrtho FCCPrim BCCPrim A/B/CCenterPrim
my $BravaisLattice = 'FCCPrim';
#my $InputMode      = 'Cartesian';
my $InputMode     = 'Reciprocal';
my $Direction;
if($InputMode eq 'Cartesian') {
	$Direction = 'OriginalToConverted';
}
else {
	$Direction = 'ConvertedlToOriginal';
}

#FCC Cartesian
my @KList;
if($InputMode eq 'Cartesian') {
	@KList = (
		[0.0,  0.0,  0.0, 'GAMMA'],
		[0.0,  0.0,  1.0, 'X'],
		[0.5,  0.0,  1.0, 'W'],
		[0.5,  0.5,  0.5, 'L'],
		[0.75, 0.75, 0.0, 'K'],
		[0.5,  0.0,  0.0, ''],
		);
}
else {
@KList = (
		[0.0,  0.0,  0.0,  'GAMMA'],
		[0.5,  0.5,  0.0,  'X'],
		[0.5,  0.75, 1.25, 'W'],
		[0.5,  0.5,  0.5,  'L'],
#		[0.375, 0.75, 0.375,  'W'],
		);
}

# 1st: No P1 RhombHex HexRhomb HexOrtho FCCPrim BCCPrim A/B/CCenterPrim
# 2nd: OriginalToConverted ConvertedlToOriginal
# 3rd: ai xyz a*i hkl
my ($sT, $T, $tRT, $RT, $tR) = Crystal->GetLatticeConversionMatrix($BravaisLattice, $Direction, 'ai', 1);
#1st argで指定された変換行列  : $sT
#2nd argで選択された変換行列  : $T
#実空間座標(x,y,z)の変換行列  : $tRT = transpose(T)
#逆空間ベクトル(a*i)の変換行列: $RT  = transpose(T)
#逆空間座標(h k l)の変換行列  : $tR  = $T
print "Conv M specified by 1st arg: sT=\n", $sT;
print "Conv M selected by 2nd arg :  T=\n", $T;
print "Conv M for (x,y,z)         : RT=\n", $RT;
print "Conv M for (a*_i)          :tRT=\n", $tRT;
print "Conv M for (h k l)         : tT=\n", $tR;

#my $GT = $T; # equal to $RT?
#$GT->transpose($GT);
#$GT = $GT->inverse();
#print "GT=\n";
#print $GT;

print "\n";
print "Conversion for [$BravaisLattice]\n";
print "Input data is for [$InputMode]\n";

for(my $i = 0 ; $i < @KList ; $i++) {
	my $pK = $KList[$i];
	print " k[$i] = ($pK->[0], $pK->[1], $pK->[2]) [$pK->[3]]\n";

	my $hkl = new Math::MatrixReal(3, 1);
	$hkl = Math::MatrixReal->new_from_cols([ [$pK->[0], $pK->[1], $pK->[2]] ]);

	my $chkl = $T * $hkl;
#	my $chkl = $RT * $hkl;
	print $chkl;
}