################################################################################
#
#  kPerl Logics/Mathematics Library Group
#  Object Library for N-dimensional Kartesian Transformation Matrix
#
#  copyright (c)
#    Karol Szafranski, 2005
#  author
#    Karol Szafranski, karol.szafranski@leibniz-fli.de
#
################################################################################
#
#  DESCRIPTION
#
# - purpose and background
#   The matrix representation of a geometrical transformation defines a
#   polynome that transforms a source coordinate object to a new coordinate
#   position. The matrix fields are defined as follows, featuring a N-dimen-
#   sional transformation matrix:
#
#             ---to--------
#              x    y   ...
#
#      | x    x2x  x2y  ...
#   from y    y2x  y2y  ...
#      |...   ...  ...  ...
#      | .    02x  02y  ...
#
# - switches may be set in two ways:
#   - library-wide, i.e. all active objects are affected. Here, the way to
#     assign switches is just an assignment to the library switch hash, like:
#       $Math::GeometNdPoint::LibGlob{switch}{-debug} = 1;
#   - object-associated: assign switches via method AddSwitch().
#   Note that lirary-wide switch assignments take effect in the moment of
#   object construction. Later changes to the library properties will not
#   take effect on existing objects.
#
# - individual description of functions at the beginning of the code blocks
#
################################################################################
#
#  OBJECT OPERATORS
#
# - data type conversion
#   bool           always TRUE
#   @{}            ( [$x2x,$x2y,...], [$y2y,$y2x,...], ... [$n2x,$n2y,...] )
#                  beware! The returned data is physically bound to the object
#
# - calculation
#   x              combine transformations $a/$b, apply transformation $b to
#                  point $a
#
#  OBJECT METHODS  for external access
#
# - housekeeping
#   new            constructor, possible initializations (via ini()):
#                  A with object arg, yielding clone of an object
#                  B first arg is number of dimensions and the following args
#                    serve to initialize the value of the object
#                    - no additional args, yielding NULL transformation
#                    - with reference to 2D ARRAY to serve for an explicit
#                      matrix assignment
#   ini            object initialization, args of new()
#   AddSwitch      assign object switches
#   clone          return copy of object
#   string         represent object in compact string format using brackets
#   dim            return/set dimension of matrix
#
# - matrix representation of common transformations
#   new_null       like method new
#   new_move       return matrix representing move by vector
#   new_rotate     return matrix representing rotation at specified Kartesian
#                  axis by radian angle
#
# - matrix calculation
#   append         add subsequent matrix transformation
#   combine        return successive combination of two matrices
#
#  OBJECT DATA STRUCTURE
#  (hash)
#
#   dim            dimension
#   mat            reference to matrix array
#                  N+1 lines, keeping N elements each
#
################################################################################
#
#  FUNCTIONS, DATA
#
#   @ISA
#   overload
#   %LibGlob
#   %_LibParam
#
# - housekeeping
#   $LibGlob{switch}
#   &new  see MainLib::DefaultObjHash.pm
#   &ini
#   &AddSwitch
#   &clone
#   &string
#   &dim
#
# - matrix representation of common transformations
#   &new_null
#   &new_move
#   &new_rotate
#
# - matrix calculation
#   &append
#   &combine
#
################################################################################

package Math::GeometNdTrform;

# includes
#use strict; use warnings;  # OK 20050623
use MainLib::DefaultObjHash;
use MainLib::Data;

# inheritance
our @ISA;
push @ISA, qw(MainLib::DefaultObjHash);

# overload
use overload (
'""'   => # make the stringify operator work like for non-overloaded object references
          sub { \%{$_[0]} },
'bool' => sub { 1 },
'@{}'  => sub { $_[0]->{mat} },
'x='   => sub { $_[0]->append($_[1]) },
'x'    => sub { &combine ($_[0], $_[1]) },
);

# package-wide constants and variables
our (%LibGlob);
my %_LibParam;


################################################################################
# housekeeping
################################################################################


$LibGlob{switch} = {};


# parametric initialization
#
# INTERFACE
# - argument 0:  class or object reference
# - argument 1:  dimension
# - argument 2*: initialization argument, either:
#                - existing object (function as a copy constructor)
#                - reference to matrix array
#                - none => null transformation
#
# - library options:
#   -debug      [STD]
#
# - return val: object reference
#
# DESCRIPTION
# - this function will be called by method new(), as a constructor sub-function.
#
sub ini {
  my ($this,$dim,$arg) = @_;

  # switches
  $this->{switch} = { %{$LibGlob{switch}} };

  # Kartesian dimensions
  if (($this->{dim}=$dim) < 2) {
    die sprintf "%s. ERROR: bad dimension: %d\n", (caller(0))[3], $dim;
  }

  # construct from submitted data structure
  if (ref ($arg)) {

    # assign from refs
    if (ref($arg) eq 'ARRAY' and int(@$arg)==3) {
      $this->{mat} = &DataClone($arg);
      $this->dim($dim);
    } elsif (ref($arg) eq __PACKAGE__) {
      $this = $arg->clone();
      $this->dim($dim);
    } else {
      die sprintf "%s. ERROR: unable to interpret argument %s\n", (caller(0))[3], $arg||"''";
    }
  }

  # assign null transformation
  else {
    $this->{mat} = [];
    for (my $i=0; $i<$dim; ++$i) {
      my $pRow = [ (0) x $dim ];
      $pRow->[$i] = 1;
      push @{$this->{mat}}, $pRow;
    }
    push @{$this->{mat}}, [ (0) x $dim ];
  }

  return $this;
}


# enter object switches
#
# INTERFACE
# - return val: success status (boolean)
#
sub AddSwitch {
  my ($this, %oopt) = @_;

  # loop over switches
  while (my($key,$val) = each(%oopt)) {
    if (0) { }

    # we don't enter -TabType into object switches
    elsif ($key eq '-TabType') { next }

    # options that we just have to enter
    else {
      if (defined $val) { $this->{switch}{$key} = $val; }
      else       { delete $this->{switch}{$key}; }
    }
  }

  # return success status
  return 1;
}


# copy object
#
# INTERFACE
# - return val: object reference
#
sub clone {
  my ($this) = @_;
  my $pCopy = &DataClone({%$this})
    or die sprintf "%s. ERROR in cloning data structure\n", (caller(0))[3];
  bless ($pCopy, ref($this));
    # this way of blessing also works for derived classes
}

sub copy { &clone(@_) }


# compact string representation of object
#
# INTERFACE
# - options:
#   -fmtfloat   sprintf format substring for float conversion, default %s
# - return val: string representation
#
sub string {
  my ($this,%opt) = @_;
  my $fmtfloat = $opt{-fmtfloat} || '%s';
  my @row;
  foreach (@{$this->{mat}}) {
    push @row, join (',', map{sprintf($fmtfloat,$_)} @$_);
  }
  sprintf ('transform(%s)', join (',', map{sprintf('[%s]',$_)} @row));
}


# return/set dimension of matrix
#
# INTERFACE
# - argument 1*: dimension (assignment mode)
# - return val:  dimension
#
sub dim {
  my ($this,$arg) = @_;
  if (defined ($arg)) {
    if ($arg < 2) {
      die sprintf "%s. ERROR: bad argument: %s\n", (caller(0))[3], $arg||"''";
    }
    my $dim = int($arg);
    # shrink matrix - move addends
    if ($arg < $this->{dim}) {
      $this->{mat}[$arg] = $this->{mat}[$this->{dim}];
      $#{$this->{mat}} = $dim;
      foreach (@{$this->{mat}}) {
        $#$_ = $dim-1;
      }
    }
    # expand matrix - move addends,
    else {
      for (my $i=$this->{dim}; $i<$dim; ++$i) {
        splice @{$this->{mat}}, $i, 0, [ (0) x $dim ];
        $this->{mat}[$i][$i] = 1;
      }
      foreach (@{$this->{mat}}) {
        while (int(@$_) < $dim) { push @$_, 0 }
      }
    }
    $this->{dim} = $dim;
  }
  return $this->{dim};
}


################################################################################
# matrix representation of common transformations
################################################################################


# matrix representation of NULL transformation
#
# INTERFACE
# - argument 1: dimension
# - return val: requested matrix
#
sub new_null {
  return $_[0]->new($_[1]);
}


# matrix representation of vector transformation
#
# INTERFACE
# - argument 1+: vector initialization arguments, including dimensions
# - return val:  requested matrix
#
sub new_move {
  require Math::GeometNdPoint;
  my ($this,@ArgVec) = @_;
  my $pVec;
  if (not $pVec = eval('Math::GeometNdPoint->new(@ArgVec)')) {
    die sprintf "%s. ERROR: unable to interpret argument as a vector: %s\n", (caller(0))[3],
      join(',',@ArgVec)||"''";
  }

  # initialize as a NULL transformation and add vector transformation
  my $pNew = $this->new($pVec->dim());
  $pNew->{mat}[$pVec->dim()] = [ $pVec->coord() ];

  return $pNew;
}


# matrix representation of stretching
#
# INTERFACE
# - argument 1: dimensions
# - argument 2: stretch factor
# - options:
#   -axis       axis to stretch along, default: global
# - return val: requested matrix
#
sub new_stretch {
  my ($this,$dim,$multiplier,%opt) = @_;
  my $pNew = $this->new($dim);

  # apply stretch factor
  if (my $a=$opt{-axis}) {
    $pNew->{mat}[$a][$a] = $multiplier;
  } else {
    for (my $a=0; $a<$dim; ++$a) {
      $pNew->{mat}[$a][$a] = $multiplier;
    }
  }
  return $pNew;
}


# matrix representation of rotation
#
# INTERFACE
# - argument 0: class or object reference
# - argument 1: Kartesian dimensions
# - argument 2: affected dimension #1
# - argument 3: affected dimension #2
# - argument 4: rotation angle (radians)
# - return val: requested matrix
#
# DESCRIPTION
# - the radian angle definition starts with 0 on the right, turning counter-
#   clockwise.
#
sub new_rotate {
  my ($this,$dim,$dimA,$dimB,$angle) = @_;

  # initialize matrix according to angle argument
  my $ValSin = sin ($angle);
  my $ValCos = cos ($angle);
  my $pNew = $this->new($dim);
  $pNew->[$dimA][$dimA] =  $ValCos;
  $pNew->[$dimA][$dimB] =  $ValSin;
  $pNew->[$dimB][$dimA] = -$ValSin;
  $pNew->[$dimB][$dimB] =  $ValCos;

  return $pNew;
}


################################################################################
# matrix calculation
################################################################################


# combine with another transformation (mutator)
#
# INTERFACE
# - return val: reference to another matrix
# - return val: object reference
#
sub append {
  my ($this, $pMat2) = @_;
  my $pComb = $this->combine($pMat2);
  $this->{mat} = $pComb->{mat};
  return $this;
}

# combine with another transformation (constructor)
#
# INTERFACE
# - return val: reference to another matrix
# - return val: joint matrix (new object)
#
sub combine {
  my ($this, $pTrf2) = @_;
  if (ref($pTrf2) ne __PACKAGE__) {
    die sprintf "%s. ERROR: unable to interpret argument: %s\n", (caller(0))[3], $pTrf2||"''";
  }

  # combine matrices
  my $pTrf1 = $this->clone();
  my $dim = $this->dim();
  for (my $s=0; $s<=$dim; ++$s) {
    for (my $t=0; $t<$dim; ++$t) {
      my $cell=0;
      for (my $m=0; $m<$dim; ++$m) {
        $cell += $pTrf2->{mat}[$m][$t]*$pTrf1->{mat}[$s][$m];
      }
      $this->{mat}[$s][$t] = $cell;
    }
  }
  for (my $t=0; $t<$dim; ++$t) {
    $this->{mat}[$dim][$t] += $pTrf2->{mat}[$dim][$t];
  }

  return $this;
}


1;
# $Id: GeometNdTrform.pm,v 1.6 2018/06/05 18:02:56 szafrans Exp $
