#!/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 MyApplication;
use Sci qw($pi);
use Crystal::CIF;
use Crystal::Crystal;
use Crystal::RietanFP;

use strict;
#use warnings;

#===============================================
# Global変数
#===============================================

#===============================================
# Applicationクラス作製
#===============================================
my $App = new MyApplication;
exit if($App->Initialize() < 0);

#===============================================
# コマンドラインオプション読み込み
#===============================================
$App->AddArgument("--TolD",                "--TolD=val: Torelance for atomic position (A) [Def 0.01 A]", "0.01");
$App->AddArgument("--xc0",                 "--xc0=val [Def 0.0]",                 "0.0");
$App->AddArgument("--yc0",                 "--yc0=val [Def 0.0]",                 "0.0");
$App->AddArgument("--zc0",                 "--zc0=val [Def 0.0]",                 "0.0");
$App->AddArgument("--ScanAll",             "--ScanAll=[0|1] [Def 0]",             "0");
$App->AddArgument("--FindSiteSymmetry",    "--FindSiteSymmetry=[0|1] [Def 0]",    "0");
$App->AddArgument("--ShowCrystalInf",      "--ShowCrystalInf=[0|1] [Def 0]",      "0");
$App->AddArgument("--ShowAllCoordination", "--ShowAllCoordination=[0|1] [Def 0]", "0");
$App->AddArgument("--RCoordMax",           "--RCoordMax=val [Def 7.0]",           "7.0");
$App->AddArgument("--ShowFoundStructure",  "--ShowFoundStructure=[0|1] [Def 0]",  "0");
$App->AddArgument("--help",                "--help: Show this help");
exit 1 if($App->ReadArgs() != 1);

my $OverlapCriteriaRadius = $App->GetGetArg('TolD');
   $OverlapCriteriaRadius = 0.01 if(!defined $OverlapCriteriaRadius);
my $xc0                   = $App->GetGetArg('xc0');
   $xc0                   = 0.0 if(!defined $xc0);
my $yc0                   = $App->GetGetArg('zc0');
   $yc0                   = 0.0 if(!defined $yc0);
my $zc0                   = $App->GetGetArg('zc0');
   $zc0                   = 0.0 if(!defined $zc0);
my $ScanAll               = $App->GetGetArg('ScanAll');
   $ScanAll               = 0 if(!defined $ScanAll);
my $FindSiteSymmetry      = $App->GetGetArg('FindSiteSymmetry');
   $FindSiteSymmetry      = 0 if(!defined $FindSiteSymmetry);
my $ShowCrystalInf        = $App->GetGetArg('ShowCrystalInf');
   $ShowCrystalInf        = 0 if(!defined $ShowCrystalInf);
my $ShowAllCoordination   = $App->GetGetArg('ShowAllCoordination');
   $ShowAllCoordination   = 0 if(!defined $ShowAllCoordination);
my $RCoordMax             = $App->GetGetArg('RCoordMax');
   $RCoordMax             = 7.0 if(!defined $RCoordMax);
my $ShowFoundStructure    = $App->GetGetArg('ShowFoundStructure');
   $ShowFoundStructure    = 0 if(!defined $ShowFoundStructure);
my $CIFFile               = $App->GetGetArg(0);
   $CIFFile               = "NaCl-converted.cif" if(!defined $CIFFile);

#CIFクラスを取得
my $cif = new CIF;
#CIFファイル名を引数に与えて読み込む
unless($cif->Read($CIFFile)) {
	print "Error: Can not read [$CIFFile]\n\n";
	exit 2;
}
print "CIF                  : ", $cif->FileName(), "\n";
print "Initial center       : ($xc0, $yc0, $zc0)\n";
print "ShowCrystalInf       : $ShowCrystalInf\n";
print "ShowAllCoordination  : $ShowAllCoordination\n";
print "ShowFoundStructure   : $ShowFoundStructure\n";
print "ScanAll              : $ScanAll\n";
print "FindSiteSymmetry     : $FindSiteSymmetry\n";
print "OverlapCriteriaRadius: $OverlapCriteriaRadius A\n";
print "RCoordMax            : $RCoordMax A\n";

my $IsPrint = 0;
&FindPossibleSPGs($OverlapCriteriaRadius, $IsPrint);

exit 0;

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

#==========================================
# サブルーチン
#==========================================
sub FindPossibleSPGs
{
	my ($OverlapCriteriaRadius, $IsPrint) = @_;

#CIFクラスの内容から、Crystalクラスを作成
#作成した時点で、Metricsの計算、イオン情報(原子量など)の読み込み、
#格子体積の計算は終わっている
	my $Crystal = $cif->GetCCrystal();

	print "\n";
	my $CrystalName = $Crystal->CrystalName();
	print "CrystalName: $CrystalName\n";
	my ($a, $b, $c, $alpha, $beta, $gamma) = $Crystal->LatticeParameters();
	print "Lattice paramters: $a $b $c  $alpha $beta $gamma\n";
#	my ($g11,$g12,$g13,$g21,$g22,$g23,$g31,$g32,$g33) = $Crystal->Metrics();
#	print "Metrics: $g11 $g22 $g33  $g12 $g23 $g31\n";
#	my $vol = $Crystal->Volume();
#	print "Volume: $vol A^3\n";
#対称要素を使い、非対称単位中の原子を展開する
#この際、Densityも計算される
	$Crystal->ExpandCoordinates();
#	my $Density = $Crystal->Density();
#	print "Density: $Density g/cm^3\n";

#Crystalクラスから、SpaceGroupクラスを取り出す
#この時点で、並進対称要素の抽出、LatticeSystemの抽出は完了している
	my $SPG = $Crystal->GetCSpaceGroup();
	my $SPGName = $SPG->SPGName();
	my $iSPG    = $SPG->iSPG();
	my $LatticeSystem = $SPG->LatticeSystem();
	print "Space Group: $SPGName ($iSPG)\n";
	print "Lattice System: $LatticeSystem\n";

	my $nTranslation = $SPG->nTranslation();
	print "nTranslation: $nTranslation\n";
	if($ShowCrystalInf) {
		for(my $i = 0 ; $i < $nTranslation ; $i++) {
			my ($x,$y,$z) = $SPG->TranslationVector($i+1);
			print "   ($x, $y, $z)\n";
		}
		print "\n";
	}

	my $nSymmetryOperation = $SPG->nSymmetryOperation();
	print "nSymmetryOperation: $nSymmetryOperation\n";
	if($ShowCrystalInf) {
		for(my $i = 0 ; $i < $nSymmetryOperation ; $i++) {
			my $symop = $SPG->SymmetryOperation($i+1);
			print "   $symop\n";
	
			my ($x1,$y1,$z1,$t1,
				$x2,$y2,$z2,$t2,
				$x3,$y3,$z3,$t3)
				= $SPG->SymmetryOperationByMatrix($i+1);
			print "      $x1 $y1 $z1  $t1\n";
			print "      $x2 $y2 $z2  $t2\n";
			print "      $x3 $y3 $z3  $t3\n";
		}
		print "\n";
	}

#Crystalクラス中の、原子の種類 AtomTypeクラスのリストをとる
	my @AtomTypeList = $Crystal->GetCAtomTypeList();
	my $nAtomType = @AtomTypeList;
	print "nAtomType: $nAtomType\n";
	if($ShowCrystalInf) {
		for(my $i = 0 ; $i < $nAtomType ; $i++) {
			my $label        = $AtomTypeList[$i]->Label();
			my $atomname     = $AtomTypeList[$i]->AtomName();
			my $charge       = $AtomTypeList[$i]->Charge();
			my $AtomicNumber = $AtomTypeList[$i]->AtomicNumber();
			my $AtomicMass   = $AtomTypeList[$i]->AtomicMass();
			my $i1 = $i+1;
			print "   #$i1: $label : $atomname [$AtomicNumber] ($charge) ($AtomicMass)\n";
		}
		print "\n";
	}
	
#Crystalクラス中の、非対称単位中の原子 AsymmetricAtomSiteクラスのリストをとる
	my @AtomSiteList = $Crystal->GetCAsymmetricAtomSiteList();
	my $nAsymmetricAtomSite = @AtomSiteList;
	print "nAsymmetricAtomSite: $nAsymmetricAtomSite\n";
	if($ShowCrystalInf or $FindSiteSymmetry) {
		for(my $i = 0 ; $i < $nAsymmetricAtomSite ; $i++) {
			my $label     = $AtomSiteList[$i]->Label();
			my $type      = $AtomSiteList[$i]->AtomName();
			my ($x,$y,$z) = $AtomSiteList[$i]->Position(1);
			my $occupancy = $AtomSiteList[$i]->Occupancy();
			my $iAtomType = $Crystal->FindiAtomType($type);

			print "   $label ($type)[#$iAtomType]: ($x, $y, $z) [$occupancy]\n";

			my ($ret, $pMatchedSyms) = &FindSiteSymmetries($Crystal, $x, $y, $z, $IsPrint) if($FindSiteSymmetry);
		}
	}

	my @ExpandedAtomSiteList = $Crystal->GetCExpandedAtomSiteList();
	my $nTotalExpandedAtomSite = $Crystal->nTotalExpandedAtomSite();
	print "nExpanededAtomSite: $nTotalExpandedAtomSite\n";
	if($ShowCrystalInf) {
		for(my $i = 0 ; $i < $nTotalExpandedAtomSite ; $i++) {
			my $label     = $ExpandedAtomSiteList[$i]->Label();
			my $type      = $ExpandedAtomSiteList[$i]->AtomName();
			my ($x,$y,$z) = $ExpandedAtomSiteList[$i]->Position(1);
			my $occupancy = $ExpandedAtomSiteList[$i]->Occupancy();
			my $iAtomType = $Crystal->FindiAtomType($type);

			print "   $label ($type)[#$iAtomType]: ($x, $y, $z) [$occupancy]\n";
		}
		print "\n";
	}

	my ($TSym, $pID) = &FindTranslations($Crystal);
	if($TSym ne '') {
		print "** Translation [$TSym] is found\n";
		if($ShowCrystalInf) {
			print "  Asymmetric sites:\n";
			my $c = 0;
			for(my $i = 0 ; $i < $nTotalExpandedAtomSite ; $i++) {
				next if($pID->[$i]);
	
				my $p         = $ExpandedAtomSiteList[$i];
				my $label     = $p->Label();
				my $type      = $p->AtomName();
				my ($x,$y,$z) = $p->Position(1);
				my $occupancy = $p->Occupancy();
				my $iAtomType = $Crystal->FindiAtomType($type);
				printf "   %3d: %4s (%2s)[#%02d]: (%8.5f, %8.5f, %8.5f) [%f]\n",
							$c, $label, $type, $iAtomType, $x, $y, $z, $occupancy;
				$c++;
			}
		}
		print "\n";
	}

	$TSym = '' if($ScanAll);

	my ($xc, $yc, $zc, $pFreedom) = &FindOrigin($Crystal, $xc0, $yc0, $zc0);
	$xc = 0 if(!$pFreedom->[0]);
	$yc = 0 if(!$pFreedom->[1]);
	$zc = 0 if(!$pFreedom->[2]);
	print "   Translation vector: (", -$xc, ", ", -$yc, ", ", -$zc, ")\n";

	for(my $i = 0 ; $i < $nTotalExpandedAtomSite ; $i++) {
		my $p              = $ExpandedAtomSiteList[$i];
		my $atomname       = $p->AtomNameOnly();
		my ($x, $y, $z)    = $p->Position();
		my ($x1, $y1, $z1) = ($x-$xc, $y-$yc, $z-$zc);
		$p->SetPosition($x1, $y1, $z1);
		printf "  Atom[%3d]: %2s: (%8.5f, %8.5f, %8.5f) to (%8.5f,%8.5f,%8.5f)\n",
					$i, $atomname, $x, $y, $z,$x1,$y1,$z1 if($ShowCrystalInf);
	}
	$Crystal->SetCAsymmetricAtomSiteList(@ExpandedAtomSiteList);
#	@ExpandedAtomSiteList = $Crystal->GetCExpandedAtomSiteList();

	my ($drive, $dir, $filename, $ext, $lastdir, $filebody) = Deps::SplitFilePath($CIFFile);
	my $TranslatedCIFFile = "$filebody-Translated.cif";
	if(!CIF->new()->CreateCIFFileFromCCrystal($Crystal, $TranslatedCIFFile, 0, 'unix')) {
		$App->print("Error: Can not create [$CIFFile].\n");
		return 0;
	}

	my ($MostProbableSPG, $MostProbableiSPG, $MostProbableiSet, $MaxnSym) = (undef, -1, -1, 0);
	for(my $iSPG = 230 ; $iSPG >= 1 ; $iSPG--) {
		for(my $iSet = 1 ; ; $iSet++) {
			my ($SPG2, $IsMatched, $mess, $pSiteId, $pAsymSites) 
					= &CheckSPG($Crystal, $TSym, $pID, $iSPG, $iSet, $OverlapCriteriaRadius, 0);
			last if(!defined $SPG2);

			if($IsMatched) {
				my $SPGName2           = $SPG2->SPGName();
				my $LatticeSystem2     = $SPG2->LatticeSystem();
				my $nTranslation       = $SPG2->nTranslation();
				my $nSymmetryOperation = $SPG2->nSymmetryOperation();
				my $nTotalSymmetry     = $nTranslation * $nSymmetryOperation;
				if($MaxnSym < $nTotalSymmetry) {
					$MaxnSym = $nTotalSymmetry;
					$MostProbableSPG  = $SPG2;
					$MostProbableiSPG = $iSPG;
					$MostProbableiSet = $iSet;
				}

				print "SPG $SPGName2 (#$iSPG-$iSet LS=$LatticeSystem2) is matched\n";
				print "   nSym=$nSymmetryOperation  nTrans=$nTranslation\n";
				my $nASSite = @$pAsymSites;
				print "   $nASSite asymmetric sites are found.\n";

				my $NewCrystal = new Crystal();
				$NewCrystal->SetCrystalName($CrystalName);
				$NewCrystal->SetLatticeParameters($a, $b, $c, $alpha, $beta, $gamma);
				$NewCrystal->SetCSpaceGroup($SPG2);

				for(my $i = 0 ; $i < @$pAsymSites ; $i++) {
#print "SiteId=$pSiteId->[$i]\n";
#					next if(defined $pSiteId->[$i]);

					my $site      = $pAsymSites->[$i];
					my $label     = $site->Label();
					my $atomname  = $site->AtomName();
					my ($x,$y,$z) = $site->Position(1);
					my $occupancy = $site->Occupancy();
					print "   $i: $atomname: ($x, $y, $z) occ=$occupancy\n" if($ShowFoundStructure);

					$NewCrystal->AddAtomSite($label, $atomname, $x, $y, $z, $occupancy);
				}

				my $NewCIFName = $CIFFile;
				$NewCIFName =~ s/\.cif//i;
				$SPGName2 =~ s/\s+//g;
				$NewCIFName .= "-$SPGName2.cif";
				print "** Save CIF file to [$NewCIFName].\n";
				CIF->new()->CreateCIFFileFromCCrystal($NewCrystal, $NewCIFName, undef, undef, $ShowFoundStructure);
			}
			else {
#				print "SPG $SPGName2 (#$iSPG-$iSet LS=$LatticeSystem2) is not matched [$mess]\n";
			}
		}
	}

	my $SPGName2           = $MostProbableSPG->SPGName();
	my $LatticeSystem2     = $MostProbableSPG->LatticeSystem();
	my $nTranslation       = $MostProbableSPG->nTranslation();
	my $nSymmetryOperation = $MostProbableSPG->nSymmetryOperation();
	my $nTotalSymmetry     = $nTranslation * $nSymmetryOperation;
	print "\n";
	print "Most matched SPG $SPGName2 (#$MostProbableiSPG-$MostProbableiSet LS=$LatticeSystem2)\n";
	print "   nSym=$nSymmetryOperation  nTrans=$nTranslation\n";

	return;
}

sub FindOrigin
{
	my ($Crystal, $xc, $yc, $zc, $IsPrint) = @_;

	my $Rmax = $RCoordMax;
	my $kRth = 0.1;

	print "\n** Find origin from initial center ($xc, $yc, $zc)\n";

	my @MArrays = &GetSymmetryOperationArraysInCartesian();

	my ($la, $lb, $lc, $alpha, $beta, $gamma) = $Crystal->LatticeParameters();
	my @AtomSiteList = $Crystal->GetCExpandedAtomSiteList();
	my $nTotalSite   = $Crystal->nTotalExpandedAtomSite();

if(0) {
	my (@site, %pSiteArray);
	for(my $i = 0 ; $i < $nTotalSite ; $i++) {
		my $site      = $AtomSiteList[$i];
		my $atomname  = $site->AtomNameOnly();
		my ($x,$y,$z) = $site->Position(1);
		my $occupancy = $site->Occupancy();
		my $ID        = sprintf("%s%04d", $atomname, int($occupancy*1000.0+1e-3));

		my %h = (
			ID          => $ID,
			'name'      => $atomname,
			'x'         => $x,
			'y'         => $y,
			'z'         => $z,
			'occupancy' => $occupancy,
			);

		$pSiteArray{$ID} = [] if(!$pSiteArray{$ID});
		my $p = $pSiteArray{$ID};
		push(@$p, \%h);
#		print "   $i: $ID ($x, $y, $z)\n";
	}

	foreach my $key (sort keys %pSiteArray) {
		my $p = $pSiteArray{$key};
		print "$key: \n";
		for(my $i = 0 ; $i < @$p ; $i++) {
			print "   ($p->[$i]{x}, $p->[$i]{y}, $p->[$i]{z}) occ=$p->[$i]{occupancy}\n";
		}
	}
}

	my %pRArray;

	my $nxmax = int($Rmax / $la) + 1;
	my $nymax = int($Rmax / $lb) + 1;
	my $nzmax = int($Rmax / $lc) + 1;
	for(my $nz = -$nzmax ; $nz <= $nzmax ; $nz++) {
	for(my $ny = -$nymax ; $ny <= $nymax ; $ny++) {
	for(my $nx = -$nxmax ; $nx <= $nxmax ; $nx++) {
		for(my $i = 0 ; $i < $nTotalSite ; $i++) {
			my $site      = $AtomSiteList[$i];
			my $atomname  = $site->AtomNameOnly();
			my ($x,$y,$z) = $site->Position(1);
			my $occupancy = $site->Occupancy();
			my $ID        = sprintf("%s%04d", $atomname, int($occupancy*1000.0+1e-3));

			my $r = $Crystal->GetInterAtomicDistance($xc, $yc, $zc, $x+$nx, $y+$ny, $z+$nz);
			next if($r > $Rmax);

			my %h = (
				ID          => $ID,
				'name'      => $atomname,
				'x'         => $x+$nx,
				'y'         => $y+$ny,
				'z'         => $z+$nz,
				'occupancy' => $occupancy,
				'r'         => $r,
				);

			$pRArray{$ID} = [] if(!$pRArray{$ID});
			my $p = $pRArray{$ID};
			push(@$p, \%h);
		}
	}
	}
	}

	my ($dx, $dy, $dz) = (0, 0, 0);
	my $c = 0;
	my @F = (1, 1, 1);
	foreach my $key (sort keys %pRArray) {
		my $p = $pRArray{$key};
		@$p = sort { $a->{r} <=> $b->{r} } @$p;

		my ($i0, $Rmin, $Rth, $n1stCoord);
		if($p->[0]{r} < 3.0*$kRth) {
			$Rmin = $p->[1]{r}-$p->[0]{r};
			$i0 = 2;
		}
		else {
			$Rmin = $p->[0]{r};
			$i0 = 1;
		}
		$Rth  = $Rmin * $kRth;
		my ($xc1, $yc1, $zc1) = (0.0, 0.0, 0.0);
		for(my $i = 0 ; $i < @$p ; $i++) {
			if($i >= $i0 and $p->[$i]{r} - $p->[$i-1]{r} > $Rth) {
				$n1stCoord = $i;
				last;
			}
#			($xc1, $yc1, $zc1) += ($p->[$i-1]{x}, $p->[$i-1]{y}, $p->[$i-1]{z});
			$xc1 += $p->[$i]{x};
			$yc1 += $p->[$i]{y};
			$zc1 += $p->[$i]{z};
		}

#		($xc1, $yc1, $zc1) /= $n1stCoord;
		if($n1stCoord > 1) {
			$xc1 /= $n1stCoord;
			$yc1 /= $n1stCoord;
			$zc1 /= $n1stCoord;
		}

		$dx += $xc1;
		$dy += $yc1;
		$dz += $zc1;
		$c++;

		printf "%7s: # of atoms in 1st coordination = %2d: Center (%8.5f, %8.5f, %8.5f)\n", 
			$key, $n1stCoord, $xc1, $yc1, $zc1;

print "    Find site symmetry at ($xc1, $yc1, $zc1)\n";
		my ($ret, $pMatchedSyms) = &FindSiteSymmetries($Crystal, $xc1, $yc1, $zc1);#, $IsPrint);
		for(my $is = 0 ; $is < @MArrays ; $is++) {
			next if(!$pMatchedSyms->[$is]);

#			my $pM   = $MArrays[$is]->{Matrix};
#			my $Name = $MArrays[$is]->{Name};
			my $pF   = $MArrays[$is]->{Freedom};
			
			$F[0] = 0 if(!$pF->[0]);
			$F[1] = 0 if(!$pF->[1]);
			$F[2] = 0 if(!$pF->[2]);
		}
print "   ** Total freedom: ($F[0], $F[1], $F[2])\n";

		if($ShowAllCoordination) {
			for(my $i = 0 ; $i < @$p ; $i++) {
				printf "   $i: (%8.5f, %8.5f, %8.5f)) r=%8.5f  occ=%f\n",
							$p->[$i]{x}, $p->[$i]{y}, $p->[$i]{z}, $p->[$i]{r}, $p->[$i]{occupancy};
			}
		}
	}
	print "\n";

	$dx /= $c;
	$dy /= $c;
	$dz /= $c;
#exit;

	return ($dx, $dy, $dz, \@F);
}

sub FindSiteSymmetries
{
	my ($Crystal, $xc, $yc, $zc, $IsPrint) = @_;

	my @MArrays = &GetSymmetryOperationArraysInCartesian();
#print "     Check site symmetry for $label [$type] at ($xc, $yc, $zc)...\n";

	my ($rxc, $ryc, $rzc) = $Crystal->CalCartesianCoordinate($xc, $yc, $zc);

	my @AtomSiteList           = $Crystal->GetCAsymmetricAtomSiteList();
	my $nAsymmetricAtomSite    = @AtomSiteList;
	my @ExpandedAtomSiteList   = $Crystal->GetCExpandedAtomSiteList();
	my $nTotalExpandedAtomSite = $Crystal->nTotalExpandedAtomSite();
#	my @nMultiplicityExpandedAtomSiteList = $Crystal->GetCnMultiplicityExpandedAtomSiteList();

	my (@hit);
	for(my $is = 0 ; $is < @MArrays ; $is++) {
		my $pM   = $MArrays[$is]->{Matrix};
		my $Name = $MArrays[$is]->{Name};
		my $pF   = $MArrays[$is]->{Freedom};

		$hit[$is] = 1;
		for(my $i = 0 ; $i < $nAsymmetricAtomSite ; $i++) {
			my $p         = $AtomSiteList[$i];
			my $label     = $p->Label();
			my $type      = $p->AtomName();
			my ($x,$y,$z) = $p->Position(1);
			my $occupancy = $p->Occupancy();

			my $x0 = $x - $xc;
			my $y0 = $y - $yc;
			my $z0 = $z - $zc;
			my ($rx0, $ry0, $rz0) = $Crystal->CalCartesianCoordinate($x0, $y0, $z0);

#			my $xs = $pM->[0][0] * $x0 + $pM->[0][1] * $y0 + $pM->[0][2] * $z0 + $pM->[0][3] + $xc;
#			my $ys = $pM->[1][0] * $x0 + $pM->[1][1] * $y0 + $pM->[1][2] * $z0 + $pM->[1][3] + $yc;
#			my $zs = $pM->[2][0] * $x0 + $pM->[2][1] * $y0 + $pM->[2][2] * $z0 + $pM->[2][3] + $zc;
			my $rxs = $pM->[0][0] * $rx0 + $pM->[0][1] * $ry0 + $pM->[0][2] * $rz0 + $pM->[0][3] + $rxc;
			my $rys = $pM->[1][0] * $rx0 + $pM->[1][1] * $ry0 + $pM->[1][2] * $rz0 + $pM->[1][3] + $ryc;
			my $rzs = $pM->[2][0] * $rx0 + $pM->[2][1] * $ry0 + $pM->[2][2] * $rz0 + $pM->[2][3] + $rzc;
			my ($xs, $ys, $zs) = $Crystal->CalFractionalCoordinate($rxs, $rys, $rzs);
#if($Name =~ /3z/) {
#print "(x,y,z)=($x,$y,$z)  0=($x0,$y0,$z0)\n";
#print "zs = $xs = $pM->[0][0] * $x0 + $pM->[0][1] * $y0 + $pM->[0][2] * $z0 + $pM->[0][3] + $xc\n";
#print "ys = $ys = $pM->[1][0] * $x0 + $pM->[1][1] * $y0 + $pM->[1][2] * $z0 + $pM->[1][3] + $yc\n";
#print "zs = $zs = $pM->[2][0] * $x0 + $pM->[2][1] * $y0 + $pM->[2][2] * $z0 + $pM->[2][3] + $zc\n";
#}

			my $id2 = $Crystal->FindIdenticalSite($type, $xs, $ys, $zs, $OverlapCriteriaRadius);#, $i);
#$IsPrint=1;
			if(defined $id2) {
printf " Sym Op[$is] ($Name): ($xs,$ys,$zs): Atom ($x,$y,$z) is equivalent to site $id2\n",
					$xs, $ys, $zs, $x, $y, $z if($IsPrint);
				next;
			}
			else {
printf " Sym Op[$is] ($Name): (%7.4f,%7.4f,%7.4f): Equivalent site of (%7.4f,%7.4f,%7.4f) was not found\n",
					$xs, $ys, $zs, $x, $y, $z if($IsPrint);
				$hit[$is] = 0;
				last;
			}
		}
	}

	my ($PG) = &FindPointGroup(\@hit, \@MArrays);
print "      Point Group: $PG\n";

print "      Matched symmetries: ";
	for(my $i = 0 ; $i < @MArrays ; $i++) {
		next if(!$hit[$i]);

		my $Name = $MArrays[$i]->{Name};
print "$Name ";
	}
print "\n";

if(0) {
print "      Incompatible symmetries: ";
	for(my $i = 0 ; $i < @MArrays ; $i++) {
		next if($hit[$i]);

		my $Name = $MArrays[$i]->{Name};
print "$Name ";
	}
print "\n";
}

	return (1, \@hit);
}

sub FindPointGroup
{
	my ($pHit, $pMArrays) = @_;

	my @MArrays;
	if(!$pMArrays) {
		@MArrays = &GetSymmetryOperationArraysInCartesian();
	}
	else {
		@MArrays = @$pMArrays;
	}
	my @PGArrays = &GetPointGroupGenerators();

	my %SymOp;
	for(my $is = 0 ; $is < @MArrays ; $is++) {
		if($pHit->[$is]) {
			$SymOp{$MArrays[$is]->{Name}} = 1;
		}
	}

	for(my $ip = 0 ; $ip < @PGArrays ; $ip++) {
		my $pGen = $PGArrays[$ip]{Generators};
		my $IsHit = 1;
		for(my $ig = 0 ; $ig < @$pGen ; $ig++) {
			if(!$SymOp{$pGen->[$ig]}) {
				$IsHit = 0;
				last;
			}
		}
		return $PGArrays[$ip]{Name} if($IsHit);
	}
	
	return '1';
}

sub FindTranslations
{
	my ($Crystal, $IsPrint) = @_;

	my @TranslationArrays = &GetTranslationArrays();

	my (@hit, $pID, @ID);
	for(my $i = 0 ; $i < @TranslationArrays ; $i++) {
		($hit[$i], $pID) = &CheckATransration($Crystal, $TranslationArrays[$i], $IsPrint);
		if($hit[$i]) {
#			print "** Translation [TranslationArrays[$i]{Name}] is found\n";
			for(my $j = 0 ; $j < @$pID ; $j++) {
				$ID[$j] = 1 if($pID->[$j]);
			}
		}
	}

	if($hit[0] and $hit[1] and $hit[2]) {
		print "** Translation [F] is found\n" if($IsPrint);
		return ("F", \@ID);
	}
	for(my $i = 0 ; $i < @TranslationArrays ; $i++) {
		if($hit[$i]) {
			print "** Translation [$TranslationArrays[$i]{Name}] is found\n" if($IsPrint);
			return ($TranslationArrays[$i]{Name}, \@ID);
		}
	}
	return ('', \@ID);
}

sub CheckATransration
{
	my ($Crystal, $ppT, $IsPrint) = @_;
#$IsPrint=1;

	my $pT = $ppT->{Matrix};

	my @SiteList = $Crystal->GetCExpandedAtomSiteList();
	my $nSite = @SiteList;

	my @SiteId;
	for(my $i = 0 ; $i < $nSite ; $i++) {
		next if($SiteId[$i]);

		my $p         = $SiteList[$i];
		my $label     = $p->Label();
		my $atomname  = $p->AtomName();
		my ($x,$y,$z) = $p->Position(1);
		my $occupancy = $p->Occupancy();
print "   $i: $label ($atomname): ($x, $y, $z) [$occupancy]\n" if($IsPrint);

		my $xs = $x + $pT->[0];
		my $ys = $y + $pT->[1];
		my $zs = $z + $pT->[2];
		my $id2 = $Crystal->FindIdenticalSite($atomname, $xs, $ys, $zs, $OverlapCriteriaRadius, $i);

		if(defined $id2) {
print " Tr Op ($pT->[0], $pT->[1], $pT->[2]): ($xs,$ys,$zs): atom site $i is equivalent to site $id2\n" if($IsPrint);
			$SiteId[$id2] = 1 if($i < $id2);
			next;
		}
		else {
print " Tr Op ($pT->[0], $pT->[1], $pT->[2]): ($xs,$ys,$zs): Equivalent site was not found for atom site $i\n" if($IsPrint);
			return (0);
		}
	}
	
	return (1, \@SiteId);
}

sub CheckSPG
{
	my ($Crystal, $TSym, $pID, $iSPG, $iSet, $OverlapCriteriaRadius, $IsPrint) = @_;

	my ($a,$b,$c,$alpha,$beta,$gamma) = $Crystal->LatticeParameters();
#	$Crystal->ExpandCoordinates();
	my $SPG = $Crystal->GetCSpaceGroup();
#	my $SPGName = $SPG->SPGName();
#	my $iSPG    = $SPG->iSPG();
	my $LatticeSystem = $SPG->LatticeSystem();
#	print "Space Group: $SPGName ($iSPG)\n";
#	print "Lattice System: $LatticeSystem\n";
	my @ExpandedAtomSiteList = $Crystal->GetCExpandedAtomSiteList();
	my $nTotalExpandedAtomSite = $Crystal->nTotalExpandedAtomSite();
	my @nMultiplicityExpandedAtomSiteList = $Crystal->GetCnMultiplicityExpandedAtomSiteList();
#	print "nTotalExpandedAtomSite: $nTotalExpandedAtomSite\n";

	my $r              = new RietanFP;
	my $SPG2           = $r->ReadSpaceGroup($iSPG, $iSet);
#print "i=$iSPG, $iSet SPG2=$SPG2\n";
	return (undef) if(!$SPG2);
	my $SPGName        = $SPG2->SPGName();
	if($TSym ne '' and $SPGName !~ /^$TSym/i) {
		return ($SPG2, 0, "[$SPGName] in SPG is not compatible with translation symmetry [$TSym]");
	}
	my $LatticeSystem2 = $SPG2->LatticeSystem();

#print "LS: [$LatticeSystem][$LatticeSystem2]\n";
	if(!&IsCompatibleLatticeSystem($LatticeSystem, $LatticeSystem2)) {
		return ($SPG2, 0, "[$LatticeSystem2] in SPG is not compatible with [$LatticeSystem]");
	}

#非対称構造中のパラメータ
	my $nTranslation       = $SPG2->nTranslation();
	my $nSymmetryOperation = $SPG2->nSymmetryOperation();
	my $nTotalSymmetry     = $nTranslation * $nSymmetryOperation;
print "nTotalSymmetry: $nTotalSymmetry = $nSymmetryOperation * $nTranslation\n" if($IsPrint);

	my @SiteId;
	for(my $i = 0 ; $i < $nTotalExpandedAtomSite ; $i++) {
		next if($pID->[$i]);
		next if($SiteId[$i]);

		my $label     = $ExpandedAtomSiteList[$i]->Label();
		my $atomname  = $ExpandedAtomSiteList[$i]->AtomName();
		my ($x,$y,$z) = $ExpandedAtomSiteList[$i]->Position(1);
		my $occupancy = $ExpandedAtomSiteList[$i]->Occupancy();
		my $id        = $ExpandedAtomSiteList[$i]->IdAsymmetricAtomSite();
		my $mult      = $nMultiplicityExpandedAtomSiteList[$id-1];
print "   $i: [$id]$label ($atomname): ($x, $y, $z) [$occupancy][$mult]\n" if($IsPrint);

		my $Found = 0;
		for(my $is = 0 ; $is < $nTotalSymmetry ; $is++) {
			my ($xs, $ys, $zs) = $SPG2->DoSymmetryOperation($is, $x, $y, $z, 1);
#print "        ($is) xs,ys,zs: $xs,$ys,$zs<br>\n";
#			my ($fxs, $fys, $fzs) = $SPG->DoVelocitySymmetryOperation($is, $fx, $fy, $fz, 1);
#			my ($vxs, $vys, $vzs) = $SPG->DoVelocitySymmetryOperation($is, $vx, $vy, $vz, 1);

#			my $id2 = $Crystal->FindIdenticalSite($atomname, $xs, $ys, $zs, $OverlapCriteriaRadius);
			my $id2 = $Crystal->FindIdenticalSite($atomname, $xs, $ys, $zs, $OverlapCriteriaRadius, $i);
			if(defined $id2) {
print " Sym Op $is ($xs,$ys,$zs): atom site $i is equivalent to site $id2\n" if($IsPrint);
				$SiteId[$id2] = 1 if($i < $id2);
				$Found = 1;
				next;
			}
			else {
print " Sym Op $is ($xs,$ys,$zs): Equivalent site was not found for atom site $i\n" if($IsPrint);
				return ($SPG2, 0, "No equivalent site found for site $i");
			}
		}
		if(!$Found) {
print "   .Equivalent site was not found for atom site $i\n" if($IsPrint);
			return ($SPG2, 0, "No equivalent site found for site $i");
		}
	}

	my @AsymSites;
	for(my $i = 0 ; $i < $nTotalExpandedAtomSite ; $i++) {
		next if(defined $SiteId[$i]);
		push(@AsymSites, $ExpandedAtomSiteList[$i]);
	}

	return ($SPG2, 1, "OK", \@SiteId, \@AsymSites);
}

sub IsCompatibleLatticeSystem
{
	my ($CrystalLS, $SPGLS) = @_;

	if($CrystalLS =~ /cub/i) {
		return 1;
	}
	if($CrystalLS =~ /tet/i) {
		return 0 if($SPGLS =~ /cub/i);
		return 1;
	}
	if($CrystalLS =~ /ort/i) {
		return 0 if($SPGLS =~ /cub/i);
		return 0 if($SPGLS =~ /tet/i);
		return 1;
	}
	if($CrystalLS =~ /hex/i) {
		return 0 if($SPGLS =~ /cub/i);
		return 0 if($SPGLS =~ /tet/i);
		return 0 if($SPGLS =~ /ort/i);
		return 1;
	}
	if($CrystalLS =~ /trig/i) {
		return 0 if($SPGLS =~ /cub/i);
		return 0 if($SPGLS =~ /tet/i);
		return 0 if($SPGLS =~ /ort/i);
		return 1;
	}
	if($CrystalLS =~ /rho/i) {
		return 0 if($SPGLS =~ /cub/i);
		return 0 if($SPGLS =~ /tet/i);
		return 0 if($SPGLS =~ /ort/i);
		return 1;
	}
	if($CrystalLS =~ /mon/i) {
		return 0 if($SPGLS =~ /cub/i);
		return 0 if($SPGLS =~ /tet/i);
		return 0 if($SPGLS =~ /ort/i);
		return 0 if($SPGLS =~ /hex/i);
		return 0 if($SPGLS =~ /trig/i);
		return 0 if($SPGLS =~ /rho/i);
		return 1;
	}
	return 1 if($SPGLS =~ /trig/i or $SPGLS eq '');
	return 0;
}

sub GetTranslationArrays
{
	my @TranslationArrays = (
		{ 'Matrix' => [0.5, 0.5, 0.0], 'Name' => 'C' },
		{ 'Matrix' => [0.5, 0.0, 0.5], 'Name' => 'B' },
		{ 'Matrix' => [0.0, 0.5, 0.5], 'Name' => 'A' },
		{ 'Matrix' => [0.5, 0.5, 0.5], 'Name' => 'I' },
		);
	
	return @TranslationArrays;
}

sub GetSymmetryOperationArraysInCartesian
{
	my $sin60  = sin( 60.0 * $pi / 180.0);
	my $cos60  = cos( 60.0 * $pi / 180.0);
	my $sin120 = sin(120.0 * $pi / 180.0);
	my $cos120 = cos(120.0 * $pi / 180.0);

	my @MArrays = (
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0, -1,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "i",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  1,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "mx",
		  'Freedom'  => [0, 1, 1],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0, -1,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "my",
		  'Freedom'  => [1, 0, 1],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0,  1,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "mz",
		  'Freedom'  => [1, 1, 0],
		},
		{ 'Matrix' => [ [ 0,  1,  0, 0],
					    [ 1,  0,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "dxy",
		  'Freedom'  => [1, 1, 1],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0,  0,  1, 0],
					    [ 0,  1,  0, 0] ],
		  'Name'     => "dyz",
		  'Freedom'  => [1, 1, 1],
		},
		{ 'Matrix' => [ [ 0,  0,  1, 0],
					    [ 0,  1,  0, 0],
					    [ 1,  0,  0, 0] ],
		  'Name'     => "dzx",
		  'Freedom'  => [1, 1, 1],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0, -1,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "2x",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  1,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "2y",
		  'Freedom'  => [0, 1, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0, -1,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "2z",
		  'Freedom'  => [0, 0, 1],
		},
#		{ 'Matrix' => [ [-1,  0,  0, 0],
#					    [ 0, -1,  0, 0],
#					    [ 0,  0, -1, 0] ],
#		  'Name'     => "S2x",
#		  'Identical' => 'i',
#		  'Freedom'  => [0, 0, 0],
#		},
#		{ 'Matrix' => [ [-1,  0,  0, 0],
#					    [ 0, -1,  0, 0],
#					    [ 0,  0, -1, 0] ],
#		  'Name'     => "S2y",
#		  'Freedom'  => [0, 0, 0],
#		},
#		{ 'Matrix' => [ [-1,  0,  0, 0],
#					    [ 0, -1,  0, 0],
#					    [ 0,  0, -1, 0] ],
#		  'Name'     => "S2z",
#		  'Identical' => 'i',
#		  'Freedom'  => [0, 0, 0],
#		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  1,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "-2x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0, -1,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "-2y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0,  1,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "-2z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [-1,  0,  1, 0],
					    [ 0,  1,  0, 0] ],
		  'Name'     => "2d011",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  0,  1, 0],
					    [ 0,  1,  0, 0] ],
		  'Name'     => "2d0-11",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0,  1, 0],
					    [ 0, -1,  0, 0],
					    [ 1,  0,  0, 0] ],
		  'Name'     => "2d101",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0, -1, 0],
					    [ 0, -1,  0, 0],
					    [-1,  0,  0, 0] ],
		  'Name'     => "2d-101",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  1,  0, 0],
					    [ 1,  0,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "2d110",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, -1,  0, 0],
					    [-1,  0,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "2d-110",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [      1,       0,       0,  0],
					    [      0, $cos120,-$sin120,  0],
					    [      0, $sin120, $cos120,  0] ],
		  'Name'     => "3x",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [ $sin120,      0, $cos120, 0],
					    [       0,      1,       0, 0],
					    [ $cos120,      0,-$sin120, 0] ],
		  'Name'     => "3y",
		  'Freedom'  => [0, 1, 0],
		},
		{ 'Matrix' => [ [ $cos120,-$sin120,      0, 0],
					    [ $sin120, $cos120,      0, 0],
					    [       0,       0,      1, 0] ],
		  'Name'     => "3z",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [      1,       0,       0,  0],
					    [      0, $cos120, $sin120,  0],
					    [      0,-$sin120, $cos120,  0] ],
		  'Name'     => "3x^2",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [-$sin120,      0, $cos120, 0],
					    [       0,      1,       0, 0],
					    [ $cos120,      0, $sin120, 0] ],
		  'Name'     => "3y^2",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [ $cos120, $sin120,      0, 0],
					    [-$sin120, $cos120,      0, 0],
					    [       0,       0,      1, 0] ],
		  'Name'     => "3z^2",
		  'Freedom'  => [0, 0, 1],
		},
		{ 'Matrix' => [ [     -1,       0,       0,  0],
					    [      0, $cos120,-$sin120,  0],
					    [      0, $sin120, $cos120,  0] ],
		  'Name'     => "S3x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $sin120,      0, $cos120, 0],
					    [       0,     -1,       0, 0],
					    [ $cos120,      0,-$sin120, 0] ],
		  'Name'     => "S3y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $cos120,-$sin120,      0, 0],
					    [ $sin120, $cos120,      0, 0],
					    [       0,       0,     -1, 0] ],
		  'Name'     => "S3z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [     -1,       0,       0,  0],
					    [      0, $cos120,-$sin120,  0],
					    [      0, $sin120, $cos120,  0] ],
		  'Name'     => "S3x^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$sin120,      0, $cos120, 0],
					    [       0,     -1,       0, 0],
					    [ $cos120,      0, $sin120, 0] ],
		  'Name'     => "S3y^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $cos120, $sin120,      0, 0],
					    [-$sin120, $cos120,      0, 0],
					    [       0,       0,     -1, 0] ],
		  'Name'     => "S3z^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [     -1,       0,       0,  0],
					    [      0,-$cos120, $sin120,  0],
					    [      0,-$sin120,-$cos120,  0] ],
		  'Name'     => "-3x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$sin120,      0,-$cos120, 0],
					    [       0,     -1,       0, 0],
					    [-$cos120,      0, $sin120, 0] ],
		  'Name'     => "-3y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$cos120, $sin120,      0, 0],
					    [-$sin120,-$cos120,      0, 0],
					    [       0,       0,     -1, 0] ],
		  'Name'     => "-3z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [     -1,       0,       0,  0],
					    [      0,-$cos120,-$sin120,  0],
					    [      0, $sin120,-$cos120,  0] ],
		  'Name'     => "-3x^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $sin120,      0,-$cos120, 0],
					    [       0,     -1,       0, 0],
					    [-$cos120,      0,-$sin120, 0] ],
		  'Name'     => "-3y^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$cos120,-$sin120,      0, 0],
					    [ $sin120,-$cos120,      0, 0],
					    [       0,       0,     -1, 0] ],
		  'Name'     => "-3z^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0,  0, -1, 0],
					    [ 0,  1,  0, 0] ],
		  'Name'     => "4x",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0,  1, 0],
					    [ 0,  1,  0, 0],
					    [-1,  0,  0, 0] ],
		  'Name'     => "4y",
		  'Freedom'  => [0, 1, 0],
		},
		{ 'Matrix' => [ [ 0, -1,  0, 0],
					    [ 1,  0,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "4z",
		  'Freedom'  => [0, 0, 1],
		},
		{ 'Matrix' => [ [ 1,  0,  0, 0],
					    [ 0,  0,  1, 0],
					    [ 0, -1,  0, 0] ],
		  'Name'     => "4x^3",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0, -1, 0],
					    [ 0,  1,  0, 0],
					    [ 1,  0,  0, 0] ],
		  'Name'     => "4y^3",
		  'Freedom'  => [0, 1, 0],
		},
		{ 'Matrix' => [ [ 0,  1,  0, 0],
					    [-1,  0,  0, 0],
					    [ 0,  0,  1, 0] ],
		  'Name'     => "4z^3",
		  'Freedom'  => [0, 0, 1],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  0, -1, 0],
					    [ 0,  1,  0, 0] ],
		  'Name'     => "S4x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0,  1, 0],
					    [ 0, -1,  0, 0],
					    [-1,  0,  0, 0] ],
		  'Name'     => "S4y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, -1,  0, 0],
					    [ 1,  0,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "S4z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  0,  1, 0],
					    [ 0, -1,  0, 0] ],
		  'Name'     => "S4x^3",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0, -1, 0],
					    [ 0, -1,  0, 0],
					    [ 1,  0,  0, 0] ],
		  'Name'     => "S4y^3",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  1,  0, 0],
					    [-1,  0,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "S4z^3",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  0,  1, 0],
					    [ 0, -1,  0, 0] ],
		  'Name'     => "-4x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0, -1, 0],
					    [ 0, -1,  0, 0],
					    [ 1,  0,  0, 0] ],
		  'Name'     => "-4y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  1,  0, 0],
					    [-1,  0,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "-4z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-1,  0,  0, 0],
					    [ 0,  0, -1, 0],
					    [ 0,  1,  0, 0] ],
		  'Name'     => "-4x^3",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0,  1, 0],
					    [ 0, -1,  0, 0],
					    [-1,  0,  0, 0] ],
		  'Name'     => "-4y^3",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, -1,  0, 0],
					    [ 1,  0,  0, 0],
					    [ 0,  0, -1, 0] ],
		  'Name'     => "-4z^3",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [      1,      0,      0,  0],
					    [      0, $cos60,-$sin60,  0],
					    [      0, $sin60, $cos60,  0] ],
		  'Name'     => "6x",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [ $sin60,      0, $cos60, 0],
					    [      0,      1,      0, 0],
					    [ $cos60,      0,-$sin60, 0] ],
		  'Name'     => "6y",
		  'Freedom'  => [0, 1, 0],
		},
		{ 'Matrix' => [ [ $cos60,-$sin60,      0, 0],
					    [ $sin60, $cos60,      0, 0],
					    [      0,      0,      1, 0] ],
		  'Name'     => "6z",
		  'Freedom'  => [0, 0, 1],
		},
		{ 'Matrix' => [ [      1,      0,      0,  0],
					    [      0, $cos60, $sin60,  0],
					    [      0,-$sin60, $cos60,  0] ],
		  'Name'     => "6x^5",
		  'Freedom'  => [1, 0, 0],
		},
		{ 'Matrix' => [ [ $sin60,      0, $cos60, 0],
					    [      0,      1,      0, 0],
					    [ $cos60,      0,-$sin60, 0] ],
		  'Name'     => "6y^5",
		  'Freedom'  => [0, 1, 0],
		},
		{ 'Matrix' => [ [ $cos60, $sin60,      0, 0],
					    [-$sin60, $cos60,      0, 0],
					    [      0,      0,      1, 0] ],
		  'Name'     => "6z^5",
		  'Freedom'  => [0, 0, 1],
		},
		{ 'Matrix' => [ [     -1,      0,      0,  0],
					    [      0, $cos60,-$sin60,  0],
					    [      0, $sin60, $cos60,  0] ],
		  'Name'     => "S6x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $sin60,      0, $cos60, 0],
					    [      0,     -1,      0, 0],
					    [ $cos60,      0,-$sin60, 0] ],
		  'Name'     => "S6y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $cos60,-$sin60,      0, 0],
					    [ $sin60, $cos60,      0, 0],
					    [      0,      0,     -1, 0] ],
		  'Name'     => "S6z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [     -1,      0,      0,  0],
					    [      0, $cos60, $sin60,  0],
					    [      0,-$sin60, $cos60,  0] ],
		  'Name'     => "S6x^5",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $sin60,      0, $cos60, 0],
					    [      0,     -1,      0, 0],
					    [ $cos60,      0,-$sin60, 0] ],
		  'Name'     => "S6y^5",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ $cos60, $sin60,      0, 0],
					    [-$sin60, $cos60,      0, 0],
					    [      0,      0,     -1, 0] ],
		  'Name'     => "S6z^5",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [     -1,      0,      0,  0],
					    [      0,-$cos60, $sin60,  0],
					    [      0,-$sin60,-$cos60,  0] ],
		  'Name'     => "-6x",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$sin60,      0,-$cos60, 0],
					    [      0,     -1,      0, 0],
					    [-$cos60,      0, $sin60, 0] ],
		  'Name'     => "-6y",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$cos60, $sin60,      0, 0],
					    [-$sin60,-$cos60,      0, 0],
					    [      0,      0,     -1, 0] ],
		  'Name'     => "-6z",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [     -1,      0,      0,  0],
					    [      0,-$cos60,-$sin60,  0],
					    [      0, $sin60,-$cos60,  0] ],
		  'Name'     => "-6x^5",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$sin60,      0,-$cos60, 0],
					    [      0,     -1,      0, 0],
					    [-$cos60,      0, $sin60, 0] ],
		  'Name'     => "-6y^5",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [-$cos60,-$sin60,      0, 0],
					    [ $sin60,-$cos60,      0, 0],
					    [      0,      0,     -1, 0] ],
		  'Name'     => "-6z^5",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, 0, 1, 0],
					    [ 1, 0, 0, 0],
					    [ 0, 1, 0, 0] ],
		  'Name'     => "3'_111",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, 1, 0, 0],
					    [ 0, 0, 1, 0],
					    [ 1, 0, 0, 0] ],
		  'Name'     => "3'111^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,-1, 0, 0],
					    [ 0, 0, 1, 0],
					    [-1, 0, 0, 0] ],
		  'Name'     => "3'_-111",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, 0,-1, 0],
					    [-1, 0, 0, 0],
					    [ 0, 1, 0, 0] ],
		  'Name'     => "3'_-111^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,-1, 0, 0],
					    [ 0, 0, 1, 0],
					    [-1, 0, 0, 0] ],
		  'Name'     => "3'_1-11",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0, 0, 1, 0],
					    [-1, 0, 0, 0],
					    [ 0,-1, 0, 0] ],
		  'Name'     => "3'_1-11^2",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  1,  0, 0],
					    [ 0,  0, -1, 0],
					    [-1,  0,  0, 0] ],
		  'Name'     => "3'_11-1",
		  'Freedom'  => [0, 0, 0],
		},
		{ 'Matrix' => [ [ 0,  0, -1, 0],
					    [ 1,  0,  0, 0],
					    [ 0, -1,  0, 0] ],
		  'Name'     => "3'_11-1^2",
		  'Freedom'  => [0, 0, 0],
		},
		);

	return @MArrays;
}

sub GetPointGroupGenerators
{
	my @PGArrays = (
		{ 'Name'       => "Oh",
		  'Generators' => ["3'_111", "3'_-111", "3'_1-11", "3'_11-1", '4x', '2x', '2d110'],
		},
		{ 'Name'       => "O",
		  'Generators' => ["3'_111", "3'_-111", "3'_1-11", "3'_11-1", '4x', '2x', '2d110', 'i', 'mx', 'dxy'],
		},
		{ 'Name'       => "Td",
		  'Generators' => ["3'_111", "3'_-111", "3'_1-11", "3'_11-1", 'S4x', 'dxy'],
		},
		{ 'Name'       => "Th",
		  'Generators' => ["3'_111", "3'_-111", "3'_1-11", "3'_11-1", '2x', 'i'],
		},
		{ 'Name'       => "T",
		  'Generators' => ["3'_111", "3'_-111", "3'_1-11", "3'_11-1", '2x'],
		},
		{ 'Name'       => "S6",
		  'Generators' => ['S6z'],
		},
		{ 'Name'       => "S6",
		  'Generators' => ['S6x'],
		},
		{ 'Name'       => "S6",
		  'Generators' => ['S6y'],
		},
		{ 'Name'       => "S4",
		  'Generators' => ['S4z'],
		},
		{ 'Name'       => "S4",
		  'Generators' => ['S4x'],
		},
		{ 'Name'       => "S4",
		  'Generators' => ['S4y'],
		},
		{ 'Name'       => "S3",
		  'Generators' => ['S3z'],
		},
		{ 'Name'       => "S3",
		  'Generators' => ['S3x'],
		},
		{ 'Name'       => "S3",
		  'Generators' => ['S3y'],
		},
		{ 'Name'       => "S2",
		  'Generators' => ['S2z'],
		},
		{ 'Name'       => "S2",
		  'Generators' => ['S2x'],
		},
		{ 'Name'       => "S2",
		  'Generators' => ['S2y'],
		},
		{ 'Name'       => "D6d",
		  'Generators' => ['6z', 'S6z', '2x', 'dxy'],
		},
		{ 'Name'       => "D6d",
		  'Generators' => ['6z', 'S6z', '2y', 'dxy'],
		},
		{ 'Name'       => "D6d",
		  'Generators' => ['6x', 'S6x', '2y', 'dyz'],
		},
		{ 'Name'       => "D6d",
		  'Generators' => ['6x', 'S6x', '2z', 'dyz'],
		},
		{ 'Name'       => "D6d",
		  'Generators' => ['6y', 'S6y', '2z', 'dzx'],
		},
		{ 'Name'       => "D6d",
		  'Generators' => ['6y', 'S6y', '2x', 'dzx'],
		},
		{ 'Name'       => "D4d",
		  'Generators' => ['4z', 'S4z', '2x', '2y', 'dxy'],
		},
		{ 'Name'       => "D4d",
		  'Generators' => ['4x', 'S4x', '2y', '2z', 'dyz'],
		},
		{ 'Name'       => "D4d",
		  'Generators' => ['4y', 'S4y', '2z', '2x', 'dzx'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3z', 'S2z', '2x', 'dxy'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3z', 'S6z', '2x', 'dxy'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3z', 'S6z', '2y', 'dxy'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3x', 'S6x', '2y', 'dyz'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3x', 'S6x', '2z', 'dyz'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3y', 'S6y', '2z', 'dzx'],
		},
		{ 'Name'       => "D3d",
		  'Generators' => ['3y', 'S6y', '2x', 'dzx'],
		},
		{ 'Name'       => "D2d",
		  'Generators' => ['2z', 'S2z', '2x', 'dxy'],
		},
		{ 'Name'       => "D2d",
		  'Generators' => ['2z', 'S2z', '2y', 'dxy'],
		},
		{ 'Name'       => "D2d",
		  'Generators' => ['2x', 'S2x', '2y', 'dyz'],
		},
		{ 'Name'       => "D2d",
		  'Generators' => ['2x', 'S2x', '2z', 'dyz'],
		},
		{ 'Name'       => "D2d",
		  'Generators' => ['2y', 'S2y', '2z', 'dzx'],
		},
		{ 'Name'       => "D2d",
		  'Generators' => ['2y', 'S2y', '2x', 'dzx'],
		},
		{ 'Name'       => "D6h",
		  'Generators' => ['6z', '2x', 'mz'],
		},
		{ 'Name'       => "D6h",
		  'Generators' => ['6z', '2y', 'mz'],
		},
		{ 'Name'       => "D6h",
		  'Generators' => ['6x', '2y', 'mx'],
		},
		{ 'Name'       => "D6h",
		  'Generators' => ['6x', '2z', 'mx'],
		},
		{ 'Name'       => "D6h",
		  'Generators' => ['6y', '2z', 'my'],
		},
		{ 'Name'       => "D6h",
		  'Generators' => ['6y', '2x', 'my'],
		},
		{ 'Name'       => "D4h",
		  'Generators' => ['4z', '2x', '2y', 'mz'],
		},
		{ 'Name'       => "D4h",
		  'Generators' => ['4x', '2y', '2z', 'mx'],
		},
		{ 'Name'       => "D4h",
		  'Generators' => ['4y', '2z', '2x', 'my'],
		},
		{ 'Name'       => "D3h",
		  'Generators' => ['3z', '2x', 'mz'],
		},
		{ 'Name'       => "D3h",
		  'Generators' => ['3z', '2y', 'mz'],
		},
		{ 'Name'       => "D3h",
		  'Generators' => ['3x', '2y', 'mx'],
		},
		{ 'Name'       => "D3h",
		  'Generators' => ['3x', '2z', 'mx'],
		},
		{ 'Name'       => "D3h",
		  'Generators' => ['3y', '2z', 'my'],
		},
		{ 'Name'       => "D3h",
		  'Generators' => ['3y', '2x', 'my'],
		},
		{ 'Name'       => "D2h",
		  'Generators' => ['2z', '2x', '2y', 'mz'],
		},
		{ 'Name'       => "D2h",
		  'Generators' => ['2x', '2y', '2z', 'mx'],
		},
		{ 'Name'       => "D2h",
		  'Generators' => ['2x', '2z', '2x', 'my'],
		},
		{ 'Name'       => "C6h",
		  'Generators' => ['6z', 'mz'],
		},
		{ 'Name'       => "C6h",
		  'Generators' => ['6x', 'mx'],
		},
		{ 'Name'       => "C6h",
		  'Generators' => ['6y', 'my'],
		},
		{ 'Name'       => "C4h",
		  'Generators' => ['4z', 'mz'],
		},
		{ 'Name'       => "C4h",
		  'Generators' => ['4x', 'mx'],
		},
		{ 'Name'       => "C4h",
		  'Generators' => ['4y', 'my'],
		},
		{ 'Name'       => "C3h",
		  'Generators' => ['3z', 'mz'],
		},
		{ 'Name'       => "C3h",
		  'Generators' => ['3x', 'mx'],
		},
		{ 'Name'       => "C3h",
		  'Generators' => ['3y', 'my'],
		},
		{ 'Name'       => "C2h",
		  'Generators' => ['2z', 'mz'],
		},
		{ 'Name'       => "C2h",
		  'Generators' => ['2x', 'mx'],
		},
		{ 'Name'       => "C2h",
		  'Generators' => ['2y', 'my'],
		},
		{ 'Name'       => "C6v",
		  'Generators' => ['6z', 'mx'],
		},
		{ 'Name'       => "C6v",
		  'Generators' => ['6z', 'my'],
		},
		{ 'Name'       => "C6v",
		  'Generators' => ['6x', 'my'],
		},
		{ 'Name'       => "C6v",
		  'Generators' => ['6x', 'mz'],
		},
		{ 'Name'       => "C6v",
		  'Generators' => ['6y', 'mz'],
		},
		{ 'Name'       => "C6v",
		  'Generators' => ['6y', 'mx'],
		},
		{ 'Name'       => "C4v",
		  'Generators' => ['4z', 'mx', 'my'],
		},
		{ 'Name'       => "C4v",
		  'Generators' => ['4x', 'my', 'mz'],
		},
		{ 'Name'       => "C4v",
		  'Generators' => ['4y', 'my', 'mz'],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ['3z', 'mx'],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ['3z', 'my'],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ['3x', 'my'],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ['3x', 'mz'],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ['3y', 'mz'],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ['3y', 'mx'],
		},
		{ 'Name'       => "C2v",
		  'Generators' => ['2z', 'mx'],
		},
		{ 'Name'       => "C2v",
		  'Generators' => ['2z', 'my'],
		},
		{ 'Name'       => "D6",
		  'Generators' => ['6z', '2x'],
		},
		{ 'Name'       => "D6",
		  'Generators' => ['6z', '2y'],
		},
		{ 'Name'       => "D6",
		  'Generators' => ['6x', '2y'],
		},
		{ 'Name'       => "D6",
		  'Generators' => ['6x', '2z'],
		},
		{ 'Name'       => "D6",
		  'Generators' => ['6y', '2z'],
		},
		{ 'Name'       => "D6",
		  'Generators' => ['6y', '2x'],
		},
		{ 'Name'       => "D4",
		  'Generators' => ['4z', '2x', '2y'],
		},
		{ 'Name'       => "D4",
		  'Generators' => ['4x', '2y', '2z'],
		},
		{ 'Name'       => "D4",
		  'Generators' => ['4y', '2y', '2z'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ['3z', '2x'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ['3z', 'C2y'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ['3x', '2y'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ['3x', '2z'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ['3y', '2z'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ['3y', '2x'],
		},
		{ 'Name'       => "D2",
		  'Generators' => ['2z', '2x', '2y'],
		},
		{ 'Name'       => "C6",
		  'Generators' => ['6x'],
		},
		{ 'Name'       => "C6",
		  'Generators' => ['6y'],
		},
		{ 'Name'       => "C6",
		  'Generators' => ['6z'],
		},
		{ 'Name'       => "C4",
		  'Generators' => ['4x'],
		},
		{ 'Name'       => "C4",
		  'Generators' => ['4y'],
		},
		{ 'Name'       => "C4",
		  'Generators' => ['4z'],
		},
		{ 'Name'       => "C3",
		  'Generators' => ['3x'],
		},
		{ 'Name'       => "C3",
		  'Generators' => ['3y'],
		},
		{ 'Name'       => "C3",
		  'Generators' => ['3z'],
		},
		{ 'Name'       => "D3",
		  'Generators' => ["3'_111", "dxy", "2d-101"],
		},
		{ 'Name'       => "D3",
		  'Generators' => ["3'_-111", "dxy", "2d-101"],
		},
		{ 'Name'       => "D3",
		  'Generators' => ["3'_1-11", "dxy", "2d-101"],
		},
		{ 'Name'       => "D3",
		  'Generators' => ["3'_11-1", "dxy", "2d-101"],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ["3'_111", "dxy"],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ["3'_-111", "dxy"],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ["3'_1-11", "dxy"],
		},
		{ 'Name'       => "C3v",
		  'Generators' => ["3'_11-1", "dxy"],
		},
		{ 'Name'       => "C3",
		  'Generators' => ["3'_111"],
		},
		{ 'Name'       => "C3",
		  'Generators' => ["3'_-111"],
		},
		{ 'Name'       => "C3",
		  'Generators' => ["3'_1-11"],
		},
		{ 'Name'       => "C3",
		  'Generators' => ["3'_11-1"],
		},
		{ 'Name'       => "C2",
		  'Generators' => ['2x'],
		},
		{ 'Name'       => "C2",
		  'Generators' => ['2y'],
		},
		{ 'Name'       => "C2",
		  'Generators' => ['2z'],
		},
		{ 'Name'       => "Ci",
		  'Generators' => ['i'],
		},
		);

	return @PGArrays;
}