A few weeks ago I needed to rotate, shear and translate a PDL image in
PGPLOT and I wanted to use the TRANSFORM parameter to get the correct
plot window.
To this end I wrote several small subroutines that make the generation
of the correct TRANSFORM vector a lot simpler.
I've attached a simple example file showing how these work. If these
are useful enough, I'd appreciate advice on how to tidy these up and
add them to PDL PGPLOT routines.
Cheers,
Matt
Matthew Kenworthy / Assistant Astronomer / Steward Observatory
933 N. Cherry Ave. / Tucson AZ 85721 / vox 520 626 6720
==============================================================
#!/usr/bin/perl -w
use PDL;
use PDL::NiceSlice;
use PDL::Graphics::PGPLOT;
dev("/xs",2,2,{WindowWidth=>8, Aspect=>1});
# BACKGROUND
#
# in PGPLOT you can remap an input image with the TRANSFORM keyword
# apart from the trivial cases, it is not easy to generate
# the required transformation vector.
#
# To address this situation, I wrote some subroutines that produce
# 3 x 3 transform matrices for translation, rotation, and scaling.
# You can then multiply these matrices together to get your
# required transform, and then a final routine reforms the 3 x 3
# matrix into the 6 element vector that you can pass to TRANSFORM
#
# transf_tx($x, $y) # make a translation matrix
# transf_rot($rotation) # make a rotation (CCW degrees) matrix
# transf_sca($scale_x, $scale_y) # make a scale (x,y) matrix
#
# mtform($m) # reform a 3x3 matrix into a 6 element TRANSFORM vector
#
# NOTE that the order of the required transforms is in REVERSE order
# i.e. if t3 x t2 x t1, then do t1 first, then t2, then t3
# EXAMPLE
#
# generate a test image of
# a star on a noisy background
$a = grandom(50,50);
$rv = rvals($a);
$a += 30*exp(-$rv*$rv/50);
# The transformation array connect the pixel index to a
# world coordinate such that:
#
# X = tr[0] + tr[1]*i + tr[2]*j
# Y = tr[3] + tr[4]*i + tr[5]*j
# a simple translation
$t1 = transf_tx(-25,-25);
# translate back
$tt1 = transf_tx(25,25);
# rotate by 30 degrees counter clockwise
$t2 = transf_rot(30);
imag $a;
label_axes "","","Original image and coordinates";
env (-30, 30, -30, 30);
imag $a, {TRANSFORM=>mtform($t1)};
label_axes "","","Translation to origin";
# the order of the transforms are in REVERSE order
# i.e. if t3 x t2 x t1, then do t1 first, then t2, then t3
env (-30, 30, -30, 30);
imag $a, {TRANSFORM=>mtform($t2 x $t1)};
label_axes "","","Rotate by 30 degrees about the centre";
# so, bring it down, rotate 30 degrees, stretch vertically by 2,
# rotate back 30 degrees
my $ang = 30;
my $sc = 2;
my $tt = transf_rot(-$ang) x transf_sca(1,$sc) x transf_rot($ang) x
transf_tx(-25, -25);
env (-50, 50, -50, 50);
imag $a, {TRANSFORM=>mtform($tt)};
label_axes "","","Rotate by 30 degrees, stretch, rotate back";
# TRANSFORMATION ROUTINES
#
# (X) = M * (i)
# (Y) (j)
# (1) (1)
#
# M = (t1 t2 t0)
# (t4 t5 t3)
# (0 0 1)
sub mtform {
use strict;
# (t1 t2 t0) I
# (t4 t5 t3) J
# (0 0 1) 1
use PDL::NiceSlice;
my ($m) = @_;
my $rem = rotate($m,1)->flat;
return($rem(0:5));
}
sub transf_tx {
use strict;
my ($x, $y) = @_;
my $tm = pdl [ [ 1, 0, $x],
[ 0, 1, $y],
[ 0, 0, 1 ] ];
return($tm);
}
sub transf_rot {
use strict;
my ($rot) = @_;
my $PI = 4 * atan2(1,1);
my $si = sin($rot * $PI / 180.);
my $co = cos($rot * $PI / 180.);
my $tm = pdl [ [ $co, -$si, 0],
[ $si, $co, 0],
[ 0, 0, 1 ] ];
return($tm);
}
sub transf_sca {
use strict;
my ($scx, $scy) = @_;
my $tm = pdl [ [ $scx, 0, 0],
[ 0, $scy, 0],
[ 0, 0, 1 ] ];
return($tm);
}
_______________________________________________
Perldl mailing list
[email protected]
http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
