1Transform(3)          User Contributed Perl Documentation         Transform(3)
2
3
4

NAME

6       PDL::Transform - Coordinate transforms, image warping, and N-D
7       functions
8

SYNOPSIS

10       use PDL::Transform;
11
12        my $t = new PDL::Transform::<type>(<opt>)
13
14        $out = $t->apply($in)  # Apply transform to some N-vectors (Transform method)
15        $out = $in->apply($t)  # Apply transform to some N-vectors (PDL method)
16
17        $im1 = $t->map($im);   # Transform image coordinates (Transform method)
18        $im1 = $im->map($t);   # Transform image coordinates (PDL method)
19
20        $t2 = $t->compose($t1);  # compose two transforms
21        $t2 = $t x $t1;          # compose two transforms (by analogy to matrix mult.)
22
23        $t3 = $t2->inverse();    # invert a transform
24        $t3 = !$t2;              # invert a transform (by analogy to logical "not")
25

DESCRIPTION

27       PDL::Transform is a convenient way to represent coordinate
28       transformations and resample images.  It embodies functions mapping R^N
29       -> R^M, both with and without inverses.  Provision exists for
30       parametrizing functions, and for composing them.  You can use this part
31       of the Transform object to keep track of arbitrary functions mapping
32       R^N -> R^M with or without inverses.
33
34       The simplest way to use a Transform object is to transform vector data
35       between coordinate systems.  The apply method accepts a PDL whose 0th
36       dimension is coordinate index (all other dimensions are threaded over)
37       and transforms the vectors into the new coordinate system.
38
39       Transform also includes image resampling, via the map method.  You
40       define a coordinate transform using a Transform object, then use it to
41       remap an image PDL.  The output is a remapped, resampled image.
42
43       You can define and compose several transformations, then apply them all
44       at once to an image.  The image is interpolated only once, when all the
45       composed transformations are applied.
46
47       In keeping with standard practice, but somewhat counterintuitively, the
48       map engine uses the inverse transform to map coordinates FROM the
49       destination dataspace (or image plane) TO the source dataspace; hence
50       PDL::Transform keeps track of both the forward and inverse transform.
51
52       For terseness and convenience, most of the constructors are exported
53       into the current package with the name "t_<transform>", so the
54       following (for example) are synonyms:
55
56         $t = new PDL::Transform::Radial();  # Long way
57
58         $t = t_radial();                    # Short way
59
60       Several math operators are overloaded, so that you can compose and
61       invert functions with expression syntax instead of method syntax (see
62       below).
63

EXAMPLE

65       Coordinate transformations and mappings are a little counterintuitive
66       at first.  Here are some examples of transforms in action:
67
68          use PDL::Transform;
69          $x = rfits('m51.fits');   # Substitute path if necessary!
70          $ts = t_linear(Scale=>3); # Scaling transform
71
72          $w = pgwin(xs);
73          $w->imag($x);
74
75          ## Grow m51 by a factor of 3; origin is at lower left.
76          $y = $ts->map($x,{pix=>1});    # pix option uses direct pixel coord system
77          $w->imag($y);
78
79          ## Shrink m51 by a factor of 3; origin still at lower left.
80          $c = $ts->unmap($x, {pix=>1});
81          $w->imag($c);
82
83          ## Grow m51 by a factor of 3; origin is at scientific origin.
84          $d = $ts->map($x,$x->hdr);    # FITS hdr template prevents autoscaling
85          $w->imag($d);
86
87          ## Shrink m51 by a factor of 3; origin is still at sci. origin.
88          $e = $ts->unmap($x,$x->hdr);
89          $w->imag($e);
90
91          ## A no-op: shrink m51 by a factor of 3, then autoscale back to size
92          $f = $ts->map($x);            # No template causes autoscaling of output
93

OPERATOR OVERLOADS

95       '!'
96          The bang is a unary inversion operator.  It binds exactly as tightly
97          as the normal bang operator.
98
99       'x'
100          By analogy to matrix multiplication, 'x' is the compose operator, so
101          these two expressions are equivalent:
102
103            $f->inverse()->compose($g)->compose($f) # long way
104            !$f x $g x $f                           # short way
105
106          Both of those expressions are equivalent to the mathematical
107          expression f^-1 o g o f, or f^-1(g(f(x))).
108
109       '**'
110          By analogy to numeric powers, you can apply an operator a positive
111          integer number of times with the ** operator:
112
113            $f->compose($f)->compose($f)  # long way
114            $f**3                         # short way
115

INTERNALS

117       Transforms are perl hashes.  Here's a list of the meaning of each key:
118
119       func
120          Ref to a subroutine that evaluates the transformed coordinates.
121          It's called with the input coordinate, and the "params" hash.  This
122          springboarding is done via explicit ref rather than by subclassing,
123          for convenience both in coding new transforms (just add the
124          appropriate sub to the module) and in adding custom transforms at
125          run-time. Note that, if possible, new "func"s should support inplace
126          operation to save memory when the data are flagged inplace.  But
127          "func" should always return its result even when flagged to compute
128          in-place.
129
130          "func" should treat the 0th dimension of its input as a dimensional
131          index (running 0..N-1 for R^N operation) and thread over all other
132          input dimensions.
133
134       inv
135          Ref to an inverse method that reverses the transformation.  It must
136          accept the same "params" hash that the forward method accepts.  This
137          key can be left undefined in cases where there is no inverse.
138
139       idim, odim
140          Number of useful dimensions for indexing on the input and output
141          sides (ie the order of the 0th dimension of the coordinates to be
142          fed in or that come out).  If this is set to 0, then as many are
143          allocated as needed.
144
145       name
146          A shorthand name for the transformation (convenient for debugging).
147          You should plan on using UNIVERAL::isa to identify classes of
148          transformation, e.g. all linear transformations should be subclasses
149          of PDL::Transform::Linear.  That makes it easier to add smarts to,
150          e.g., the compose() method.
151
152       itype
153          An array containing the name of the quantity that is expected from
154          the input piddle for the transform, for each dimension.  This field
155          is advisory, and can be left blank if there's no obvious quantity
156          associated with the transform.  This is analogous to the CTYPEn
157          field used in FITS headers.
158
159       oname
160          Same as itype, but reporting what quantity is delivered for each
161          dimension.
162
163       iunit
164          The units expected on input, if a specific unit (e.g. degrees) is
165          expected.  This field is advisory, and can be left blank if there's
166          no obvious unit associated with the transform.
167
168       ounit
169          Same as iunit, but reporting what quantity is delivered for each
170          dimension.
171
172       params
173          Hash ref containing relevant parameters or anything else the func
174          needs to work right.
175
176       is_inverse
177          Bit indicating whether the transform has been inverted.  That is
178          useful for some stringifications (see the PDL::Transform::Linear
179          stringifier), and may be useful for other things.
180
181       Transforms should be inplace-aware where possible, to prevent excessive
182       memory usage.
183
184       If you define a new type of transform, consider generating a new
185       stringify method for it.  Just define the sub "stringify" in the
186       subclass package.  It should call SUPER::stringify to generate the
187       first line (though the PDL::Transform::Composition bends this rule by
188       tweaking the top-level line), then output (indented) additional lines
189       as necessary to fully describe the transformation.
190

NOTES

192       Transforms have a mechanism for labeling the units and type of each
193       coordinate, but it is just advisory.  A routine to identify and, if
194       necessary, modify units by scaling would be a good idea.  Currently, it
195       just assumes that the coordinates are correct for (e.g.)  FITS
196       scientific-to-pixel transformations.
197
198       Composition works OK but should probably be done in a more
199       sophisticated way so that, for example, linear transformations are
200       combined at the matrix level instead of just strung together pixel-to-
201       pixel.
202

MODULE INTERFACE

204       There are both operators and constructors.  The constructors are all
205       exported, all begin with "t_", and all return objects that are
206       subclasses of PDL::Transform.
207
208       The apply, invert, map, and unmap methods are also exported to the
209       "PDL" package: they are both Transform methods and PDL methods.
210

FUNCTIONS

212   apply
213         Signature: (data(); PDL::Transform t)
214
215         $out = $data->apply($t);
216         $out = $t->apply($data);
217
218       Apply a transformation to some input coordinates.
219
220       In the example, $t is a PDL::Transform and $data is a PDL to be
221       interpreted as a collection of N-vectors (with index in the 0th
222       dimension).  The output is a similar but transformed PDL.
223
224       For convenience, this is both a PDL method and a Transform method.
225
226   invert
227         Signature: (data(); PDL::Transform t)
228
229         $out = $t->invert($data);
230         $out = $data->invert($t);
231
232       Apply an inverse transformation to some input coordinates.
233
234       In the example, $t is a PDL::Transform and $data is a piddle to be
235       interpreted as a collection of N-vectors (with index in the 0th
236       dimension).  The output is a similar piddle.
237
238       For convenience this is both a PDL method and a PDL::Transform method.
239
240   map
241         Signature: (k0(); SV *in; SV *out; SV *map; SV *boundary; SV *method;
242                           SV *big; SV *blur; SV *sv_min; SV *flux; SV *bv)
243
244   match
245         $y = $x->match($c);                  # Match $c's header and size
246         $y = $x->match([100,200]);           # Rescale to 100x200 pixels
247         $y = $x->match([100,200],{rect=>1}); # Rescale and remove rotation/skew.
248
249       Resample a scientific image to the same coordinate system as another.
250
251       The example above is syntactic sugar for
252
253        $y = $x->map(t_identity, $c, ...);
254
255       it resamples the input PDL with the identity transformation in
256       scientific coordinates, and matches the pixel coordinate system to $c's
257       FITS header.
258
259       There is one difference between match and map: match makes the
260       "rectify" option to "map" default to 0, not 1.  This only affects
261       matching where autoscaling is required (i.e. the array ref example
262       above).  By default, that example simply scales $x to the new size and
263       maintains any rotation or skew in its scientific-to-pixel coordinate
264       transform.
265
266   map
267         $y = $x->map($xform,[<template>],[\%opt]); # Distort $x with transform $xform
268         $y = $x->map(t_identity,[$pdl],[\%opt]); # rescale $x to match $pdl's dims.
269
270       Resample an image or N-D dataset using a coordinate transform.
271
272       The data are resampled so that the new pixel indices are proportional
273       to the transformed coordinates rather than the original ones.
274
275       The operation uses the inverse transform: each output pixel location is
276       inverse-transformed back to a location in the original dataset, and the
277       value is interpolated or sampled appropriately and copied into the
278       output domain.  A variety of sampling options are available, trading
279       off speed and mathematical correctness.
280
281       For convenience, this is both a PDL method and a PDL::Transform method.
282
283       "map" is FITS-aware: if there is a FITS header in the input data, then
284       the coordinate transform acts on the scientific coordinate system
285       rather than the pixel coordinate system.
286
287       By default, the output coordinates are separated from pixel coordinates
288       by a single layer of indirection.  You can specify the mapping between
289       output transform (scientific) coordinates to pixel coordinates using
290       the "orange" and "irange" options (see below), or by supplying a FITS
291       header in the template.
292
293       If you don't specify an output transform, then the output is
294       autoscaled: "map" transforms a few vectors in the forward direction to
295       generate a mapping that will put most of the data on the image plane,
296       for most transformations.  The calculated mapping gets stuck in the
297       output's FITS header.
298
299       Autoscaling is especially useful for rescaling images -- if you specify
300       the identity transform and allow autoscaling, you duplicate the
301       functionality of rescale2d, but with more options for interpolation.
302
303       You can operate in pixel space, and avoid autoscaling of the output, by
304       setting the "nofits" option (see below).
305
306       The output has the same data type as the input.  This is a feature, but
307       it can lead to strange-looking banding behaviors if you use
308       interpolation on an integer input variable.
309
310       The "template" can be one of:
311
312       ·  a PDL
313
314          The PDL and its header are copied to the output array, which is then
315          populated with data.  If the PDL has a FITS header, then the FITS
316          transform is automatically applied so that $t applies to the output
317          scientific coordinates and not to the output pixel coordinates.  In
318          this case the NAXIS fields of the FITS header are ignored.
319
320       ·  a FITS header stored as a hash ref
321
322          The FITS NAXIS fields are used to define the output array, and the
323          FITS transformation is applied to the coordinates so that $t applies
324          to the output scientific coordinates.
325
326       ·  a list ref
327
328          This is a list of dimensions for the output array.  The code
329          estimates appropriate pixel scaling factors to fill the available
330          space.  The scaling factors are placed in the output FITS header.
331
332       ·  nothing
333
334          In this case, the input image size is used as a template, and
335          scaling is done as with the list ref case (above).
336
337       OPTIONS:
338
339       The following options are interpreted:
340
341       b, bound, boundary, Boundary (default = 'truncate')
342          This is the boundary condition to be applied to the input image; it
343          is passed verbatim to range or interpND in the sampling or
344          interpolating stage.  Other values are 'forbid','extend', and
345          'periodic'.  You can abbreviate this to a single letter.  The
346          default 'truncate' causes the entire notional space outside the
347          original image to be filled with 0.
348
349       pix, Pixel, nf, nofits, NoFITS (default = 0)
350          If you set this to a true value, then FITS headers and
351          interpretation are ignored; the transformation is treated as being
352          in raw pixel coordinates.
353
354       j, J, just, justify, Justify (default = 0)
355          If you set this to 1, then output pixels are autoscaled to have unit
356          aspect ratio in the output coordinates.  If you set it to a non-1
357          value, then it is the aspect ratio between the first dimension and
358          all subsequent dimensions -- or, for a 2-D transformation, the
359          scientific pixel aspect ratio.  Values less than 1 shrink the scale
360          in the first dimension compared to the other dimensions; values
361          greater than 1 enlarge it compared to the other dimensions.  (This
362          is the same sense as in the PGPLOT interface.)
363
364       ir, irange, input_range, Input_Range
365          This is a way to modify the autoscaling.  It specifies the range of
366          input scientific (not necessarily pixel) coordinates that you want
367          to be mapped to the output image.  It can be either a nested array
368          ref or a piddle.  The 0th dim (outside coordinate in the array ref)
369          is dimension index in the data; the 1st dim should have order 2.
370          For example, passing in either [[-1,2],[3,4]] or pdl([[-1,2],[3,4]])
371          limites the map to the quadrilateral in input space defined by the
372          four points (-1,3), (-1,4), (2,4), and (2,3).
373
374          As with plain autoscaling, the quadrilateral gets sparsely sampled
375          by the autoranger, so pathological transformations can give you
376          strange results.
377
378          This parameter is overridden by "orange", below.
379
380       or, orange, output_range, Output_Range
381          This sets the window of output space that is to be sampled onto the
382          output array.  It works exactly like "irange", except that it
383          specifies a quadrilateral in output space.  Since the output pixel
384          array is itself a quadrilateral, you get pretty much exactly what
385          you asked for.
386
387          This parameter overrides "irange", if both are specified.  It forces
388          rectification of the output (so that scientific axes align with the
389          pixel grid).
390
391       r, rect, rectify
392          This option defaults TRUE and controls how autoscaling is performed.
393          If it is true or undefined, then autoscaling adjusts so that pixel
394          coordinates in the output image are proportional to individual
395          scientific coordinates.  If it is false, then autoscaling
396          automatically applies the inverse of any input FITS transformation
397          *before* autoscaling the pixels.  In the special case of linear
398          transformations, this preserves the rectangular shape of the
399          original pixel grid and makes output pixel coordinate proportional
400          to input coordinate.
401
402       m, method, Method
403          This option controls the interpolation method to be used.
404          Interpolation greatly affects both speed and quality of output.  For
405          most cases the option is directly passed to interpND for
406          interpolation.  Possible options, in order from fastest to slowest,
407          are:
408
409          ·  s, sample (default for ints)
410
411             Pixel values in the output plane are sampled from the closest
412             data value in the input plane.  This is very fast but not very
413             accurate for either magnification or decimation (shrinking).  It
414             is the default for templates of integer type.
415
416          ·  l, linear (default for floats)
417
418             Pixel values are linearly interpolated from the closest data
419             value in the input plane.  This is reasonably fast but only
420             accurate for magnification.  Decimation (shrinking) of the image
421             causes aliasing and loss of photometry as features fall between
422             the samples.  It is the default for floating-point templates.
423
424          ·  c, cubic
425
426             Pixel values are interpolated using an N-cubic scheme from a
427             4-pixel N-cube around each coordinate value.  As with linear
428             interpolation, this is only accurate for magnification.
429
430          ·  f, fft
431
432             Pixel values are interpolated using the term coefficients of the
433             Fourier transform of the original data.  This is the most
434             appropriate technique for some kinds of data, but can yield
435             undesired "ringing" for expansion of normal images.  Best suited
436             to studying images with repetitive or wavelike features.
437
438          ·  h, hanning
439
440             Pixel values are filtered through a spatially-variable filter
441             tuned to the computed Jacobian of the transformation, with
442             hanning-window (cosine) pixel rolloff in each dimension.  This
443             prevents aliasing in the case where the image is distorted or
444             shrunk, but allows small amounts of aliasing at pixel edges
445             wherever the image is enlarged.
446
447          ·  g, gaussian, j, jacobian
448
449             Pixel values are filtered through a spatially-variable filter
450             tuned to the computed Jacobian of the transformation, with radial
451             Gaussian rolloff.  This is the most accurate resampling method,
452             in the sense of introducing the fewest artifacts into a properly
453             sampled data set.  This method uses a lookup table to speed up
454             calculation of the Gaussian weighting.
455
456          ·  G
457
458             This method works exactly like 'g' (above), except that the
459             Gaussian values are computed explicitly for every sample -- which
460             is considerably slower.
461
462       blur, Blur (default = 1.0)
463          This value scales the input-space footprint of each output pixel in
464          the gaussian and hanning methods. It's retained for historical
465          reasons.  Larger values yield blurrier images; values significantly
466          smaller than unity cause aliasing.  The parameter has slightly
467          different meanings for Gaussian and Hanning interpolation.  For
468          Hanning interpolation, numbers smaller than unity control the
469          sharpness of the border at the edge of each pixel (so that blur=>0
470          is equivalent to sampling for non-decimating transforms).  For
471          Gaussian interpolation, the blur factor parameter controls the width
472          parameter of the Gaussian around the center of each pixel.
473
474       sv, SV (default = 1.0)
475          This value lets you set the lower limit of the transformation's
476          singular values in the hanning and gaussian methods, limiting the
477          minimum radius of influence associated with each output pixel.
478          Large numbers yield smoother interpolation in magnified parts of the
479          image but don't affect reduced parts of the image.
480
481       big, Big (default = 0.5)
482          This is the largest allowable input spot size which may be mapped to
483          a single output pixel by the hanning and gaussian methods, in units
484          of the largest non-thread input dimension.  (i.e. the default won't
485          let you reduce the original image to less than 5 pixels across).
486          This places a limit on how long the processing can take for
487          pathological transformations.  Smaller numbers keep the code from
488          hanging for a long time; larger numbers provide for photometric
489          accuracy in more pathological cases.  Numbers larer than 1.0 are
490          silly, because they allow the entire input array to be compressed
491          into a region smaller than a single pixel.
492
493          Wherever an output pixel would require averaging over an area that
494          is too big in input space, it instead gets NaN or the bad value.
495
496       phot, photometry, Photometry
497          This lets you set the style of photometric conversion to be used in
498          the hanning or gaussian methods.  You may choose:
499
500          ·  0, s, surf, surface, Surface (default)
501
502             (this is the default): surface brightness is preserved over the
503             transformation, so features maintain their original intensity.
504             This is what the sampling and interpolation methods do.
505
506          ·  1, f, flux, Flux
507
508             Total flux is preserved over the transformation, so that the
509             brightness integral over image regions is preserved.  Parts of
510             the image that are shrunk wind up brighter; parts that are
511             enlarged end up fainter.
512
513       VARIABLE FILTERING:
514
515       The 'hanning' and 'gaussian' methods of interpolation give
516       photometrically accurate resampling of the input data for arbitrary
517       transformations.  At each pixel, the code generates a linear
518       approximation to the input transformation, and uses that linearization
519       to estimate the "footprint" of the output pixel in the input space.
520       The output value is a weighted average of the appropriate input spaces.
521
522       A caveat about these methods is that they assume the transformation is
523       continuous.  Transformations that contain discontinuities will give
524       incorrect results near the discontinuity.  In particular, the 180th
525       meridian isn't handled well in lat/lon mapping transformations (see
526       PDL::Transform::Cartography) -- pixels along the 180th meridian get the
527       average value of everything along the parallel occupied by the pixel.
528       This flaw is inherent in the assumptions that underly creating a
529       Jacobian matrix.  Maybe someone will write code to work around it.
530       Maybe that someone is you.
531
532       BAD VALUES:
533
534       If your PDL was compiled with bad value support, "map()" supports bad
535       values in the data array.  Bad values in the input array are propagated
536       to the output array.  The 'g' and 'h' methods will do some smoothing
537       over bad values:  if more than 1/3 of the weighted input-array
538       footprint of a given output pixel is bad, then the output pixel gets
539       marked bad.
540
541       map does not process bad values.  It will set the bad-value flag of all
542       output piddles if the flag is set for any of the input piddles.
543
544   unmap
545        Signature: (data(); PDL::Transform a; template(); \%opt)
546
547         $out_image = $in_image->unmap($t,[<options>],[<template>]);
548         $out_image = $t->unmap($in_image,[<options>],[<template>]);
549
550       Map an image or N-D dataset using the inverse as a coordinate
551       transform.
552
553       This convenience function just inverts $t and calls map on the inverse;
554       everything works the same otherwise.  For convenience, it is both a PDL
555       method and a PDL::Transform method.
556
557   t_inverse
558         $t2 = t_inverse($t);
559         $t2 = $t->inverse;
560         $t2 = $t ** -1;
561         $t2 = !$t;
562
563       Return the inverse of a PDL::Transform.  This just reverses the
564       func/inv, idim/odim, itype/otype, and iunit/ounit pairs.  Note that
565       sometimes you end up with a transform that cannot be applied or mapped,
566       because either the mathematical inverse doesn't exist or the inverse
567       func isn't implemented.
568
569       You can invert a transform by raising it to a negative power, or by
570       negating it with '!'.
571
572       The inverse transform remains connected to the main transform because
573       they both point to the original parameters hash.  That turns out to be
574       useful.
575
576   t_compose
577         $f2 = t_compose($f, $g,[...]);
578         $f2 = $f->compose($g[,$h,$i,...]);
579         $f2 = $f x $g x ...;
580
581       Function composition: f(g(x)), f(g(h(x))), ...
582
583       You can also compose transforms using the overloaded matrix-
584       multiplication (nee repeat) operator 'x'.
585
586       This is accomplished by inserting a splicing code ref into the "func"
587       and "inv" slots.  It combines multiple compositions into a single list
588       of transforms to be executed in order, fram last to first (in keeping
589       with standard mathematical notation).  If one of the functions is
590       itself a composition, it is interpolated into the list rather than left
591       separate.  Ultimately, linear transformations may also be combined
592       within the list.
593
594       No checking is done that the itype/otype and iunit/ounit fields are
595       compatible -- that may happen later, or you can implement it yourself
596       if you like.
597
598   t_wrap
599         $g1fg = $f->wrap($g);
600         $g1fg = t_wrap($f,$g);
601
602       Shift a transform into a different space by 'wrapping' it with a
603       second.
604
605       This is just a convenience function for two t_compose calls.
606       "$x->wrap($y)" is the same as "(!$y) x $x x $y": the resulting
607       transform first hits the data with $y, then with $x, then with the
608       inverse of $y.
609
610       For example, to shift the origin of rotation, do this:
611
612         $im = rfits('m51.fits');
613         $tf = t_fits($im);
614         $tr = t_linear({rot=>30});
615         $im1 = $tr->map($tr);               # Rotate around pixel origin
616         $im2 = $tr->map($tr->wrap($tf));    # Rotate round FITS scientific origin
617
618   t_identity
619         my $xform = t_identity
620         my $xform = new PDL::Transform;
621
622       Generic constructor generates the identity transform.
623
624       This constructor really is trivial -- it is mainly used by the other
625       transform constructors.  It takes no parameters and returns the
626       identity transform.
627
628   t_lookup
629         $f = t_lookup($lookup, {<options>});
630
631       Transform by lookup into an explicit table.
632
633       You specify an N+1-D PDL that is interpreted as an N-D lookup table of
634       column vectors (vector index comes last).  The last dimension has order
635       equal to the output dimensionality of the transform.
636
637       For added flexibility in data space, You can specify pre-lookup linear
638       scaling and offset of the data.  Of course you can specify the
639       interpolation method to be used.  The linear scaling stuff is a little
640       primitive; if you want more, try composing the linear transform with
641       this one.
642
643       The prescribed values in the lookup table are treated as pixel-
644       centered: that is, if your input array has N elements per row then
645       valid data exist between the locations (-0.5) and (N-0.5) in lookup
646       pixel space, because the pixels (which are numbered from 0 to N-1) are
647       centered on their locations.
648
649       Lookup is done using interpND, so the boundary conditions and threading
650       behaviour follow from that.
651
652       The indexed-over dimensions come first in the table, followed by a
653       single dimension containing the column vector to be output for each set
654       of other dimensions -- ie to output 2-vectors from 2 input parameters,
655       each of which can range from 0 to 49, you want an index that has
656       dimension list (50,50,2).  For the identity lookup table you could use
657       "cat(xvals(50,50),yvals(50,50))".
658
659       If you want to output a single value per input vector, you still need
660       that last index threading dimension -- if necessary, use "dummy(-1,1)".
661
662       The lookup index scaling is: out = lookup[ (scale * data) + offset ].
663
664       A simplistic table inversion routine is included.  This means that you
665       can (for example) use the "map" method with "t_lookup" transformations.
666       But the table inversion is exceedingly slow, and not practical for
667       tables larger than about 100x100.  The inversion table is calculated in
668       its entirety the first time it is needed, and then cached until the
669       object is destroyed.
670
671       Options are listed below; there are several synonyms for each.
672
673       s, scale, Scale
674          (default 1.0) Specifies the linear amount of scaling to be done
675          before lookup.  You can feed in a scalar or an N-vector; other
676          values may cause trouble.  If you want to save space in your table,
677          then specify smaller scale numbers.
678
679       o, offset, Offset
680          (default 0.0) Specifies the linear amount of offset before lookup.
681          This is only a scalar, because it is intended to let you switch to
682          corner-centered coordinates if you want to (just feed in o=-0.25).
683
684       b, bound, boundary, Boundary
685          Boundary condition to be fed to interpND
686
687       m, method, Method
688          Interpolation method to be fed to interpND
689
690       EXAMPLE
691
692       To scale logarithmically the Y axis of m51, try:
693
694         $x = float rfits('m51.fits');
695         $lookup = xvals(128,128) -> cat( 10**(yvals(128,128)/50) * 256/10**2.55 );
696         $t = t_lookup($lookup);
697         $y = $t->map($x);
698
699       To do the same thing but with a smaller lookup table, try:
700
701         $lookup = 16 * xvals(17,17)->cat(10**(yvals(17,17)/(100/16)) * 16/10**2.55);
702         $t = t_lookup($lookup,{scale=>1/16.0});
703         $y = $t->map($x);
704
705       (Notice that, although the lookup table coordinates are is divided by
706       16, it is a 17x17 -- so linear interpolation works right to the edge of
707       the original domain.)
708
709       NOTES
710
711       Inverses are not yet implemented -- the best way to do it might be by
712       judicious use of map() on the forward transformation.
713
714       the type/unit fields are ignored.
715
716   t_linear
717       $f = t_linear({options});
718
719       Linear (affine) transformations with optional offset
720
721       t_linear implements simple matrix multiplication with offset, also
722       known as the affine transformations.
723
724       You specify the linear transformation with pre-offset, a mixing matrix,
725       and a post-offset.  That overspecifies the transformation, so you can
726       choose your favorite method to specify the transform you want.  The
727       inverse transform is automagically generated, provided that it actually
728       exists (the transform matrix is invertible).  Otherwise, the inverse
729       transform just croaks.
730
731       Extra dimensions in the input vector are ignored, so if you pass a 3xN
732       vector into a 3-D linear transformation, the final dimension is passed
733       through unchanged.
734
735       The options you can usefully pass in are:
736
737       s, scale, Scale
738          A scaling scalar (heh), vector, or matrix.  If you specify a vector
739          it is treated as a diagonal matrix (for convenience).  It gets left-
740          multiplied with the transformation matrix you specify (or the
741          identity), so that if you specify both a scale and a matrix the
742          scaling is done after the rotation or skewing or whatever.
743
744       r, rot, rota, rotation, Rotation
745          A rotation angle in degrees -- useful for 2-D and 3-D data only.  If
746          you pass in a scalar, it specifies a rotation from the 0th axis
747          toward the 1st axis.  If you pass in a 3-vector as either a PDL or
748          an array ref (as in "rot=>[3,4,5]"), then it is treated as a set of
749          Euler angles in three dimensions, and a rotation matrix is generated
750          that does the following, in order:
751
752          ·  Rotate by rot->(2) degrees from 0th to 1st axis
753
754          ·  Rotate by rot->(1) degrees from the 2nd to the 0th axis
755
756          ·  Rotate by rot->(0) degrees from the 1st to the 2nd axis
757
758          The rotation matrix is left-multiplied with the transformation
759          matrix you specify, so that if you specify both rotation and a
760          general matrix the rotation happens after the more general operation
761          -- though that is deprecated.
762
763          Of course, you can duplicate this functionality -- and get more
764          general -- by generating your own rotation matrix and feeding it in
765          with the "matrix" option.
766
767       m, matrix, Matrix
768          The transformation matrix.  It does not even have to be square, if
769          you want to change the dimensionality of your input.  If it is
770          invertible (note: must be square for that), then you automagically
771          get an inverse transform too.
772
773       pre, preoffset, offset, Offset
774          The vector to be added to the data before they get multiplied by the
775          matrix (equivalent of CRVAL in FITS, if you are converting from
776          scientific to pixel units).
777
778       post, postoffset, shift, Shift
779          The vector to be added to the data after it gets multiplied by the
780          matrix (equivalent of CRPIX-1 in FITS, if youre converting from
781          scientific to pixel units).
782
783       d, dim, dims, Dims
784          Most of the time it is obvious how many dimensions you want to deal
785          with: if you supply a matrix, it defines the transformation; if you
786          input offset vectors in the "pre" and "post" options, those define
787          the number of dimensions.  But if you only supply scalars, there is
788          no way to tell and the default number of dimensions is 2.  This
789          provides a way to do, e.g., 3-D scaling: just set
790          "{s="<scale-factor>, dims=>3}> and you are on your way.
791
792       NOTES
793
794       the type/unit fields are currently ignored by t_linear.
795
796   t_scale
797         $f = t_scale(<scale>)
798
799       Convenience interface to t_linear.
800
801       t_scale produces a transform that scales around the origin by a fixed
802       amount.  It acts exactly the same as "t_linear(Scale="\<scale\>)>.
803
804   t_offset
805         $f = t_offset(<shift>)
806
807       Convenience interface to t_linear.
808
809       t_offset produces a transform that shifts the origin to a new location.
810       It acts exactly the same as "t_linear(Pre="\<shift\>)>.
811
812   t_rot
813         $f = t_rot(<rotation-in-degrees>)
814
815       Convenience interface to t_linear.
816
817       t_rot produces a rotation transform in 2-D (scalar), 3-D (3-vector), or
818       N-D (matrix).  It acts exactly the same as "t_linear(Rot="\<shift\>)>.
819
820   t_fits
821         $f = t_fits($fits,[option]);
822
823       FITS pixel-to-scientific transformation with inverse
824
825       You feed in a hash ref or a PDL with one of those as a header, and you
826       get back a transform that converts 0-originated, pixel-centered
827       coordinates into scientific coordinates via the transformation in the
828       FITS header.  For most FITS headers, the transform is reversible, so
829       applying the inverse goes the other way.  This is just a convenience
830       subclass of PDL::Transform::Linear, but with unit/type support using
831       the FITS header you supply.
832
833       For now, this transform is rather limited -- it really ought to accept
834       units differences and stuff like that, but they are just ignored for
835       now.  Probably that would require putting units into the whole
836       transform framework.
837
838       This transform implements the linear transform part of the WCS FITS
839       standard outlined in Greisen & Calabata 2002 (A&A in press; find it at
840       "http://arxiv.org/abs/astro-ph/0207407").
841
842       As a special case, you can pass in the boolean option "ignore_rgb"
843       (default 0), and if you pass in a 3-D FITS header in which the last
844       dimension has exactly 3 elements, it will be ignored in the output
845       transformation.  That turns out to be handy for handling rgb images.
846
847   t_code
848         $f = t_code(<func>,[<inv>],[options]);
849
850       Transform implementing arbitrary perl code.
851
852       This is a way of getting quick-and-dirty new transforms.  You pass in
853       anonymous (or otherwise) code refs pointing to subroutines that
854       implement the forward and, optionally, inverse transforms.  The
855       subroutines should accept a data PDL followed by a parameter hash ref,
856       and return the transformed data PDL.  The parameter hash ref can be set
857       via the options, if you want to.
858
859       Options that are accepted are:
860
861       p,params
862          The parameter hash that will be passed back to your code (defaults
863          to the empty hash).
864
865       n,name
866          The name of the transform (defaults to "code").
867
868       i, idim (default 2)
869          The number of input dimensions (additional ones should be passed
870          through unchanged)
871
872       o, odim (default 2)
873          The number of output dimensions
874
875       itype
876          The type of the input dimensions, in an array ref (optional and
877          advisiory)
878
879       otype
880          The type of the output dimension, in an array ref (optional and
881          advisory)
882
883       iunit
884          The units that are expected for the input dimensions (optional and
885          advisory)
886
887       ounit
888          The units that are returned in the output (optional and advisory).
889
890       The code variables are executable perl code, either as a code ref or as
891       a string that will be eval'ed to produce code refs.  If you pass in a
892       string, it gets eval'ed at call time to get a code ref.  If it compiles
893       OK but does not return a code ref, then it gets re-evaluated with "sub
894       { ... }" wrapped around it, to get a code ref.
895
896       Note that code callbacks like this can be used to do really weird
897       things and generate equally weird results -- caveat scriptor!
898
899   t_cylindrical
900       "t_cylindrical" is an alias for "t_radial"
901
902   t_radial
903         $f = t_radial(<options>);
904
905       Convert Cartesian to radial/cylindrical coordinates.  (2-D/3-D; with
906       inverse)
907
908       Converts 2-D Cartesian to radial (theta,r) coordinates.  You can choose
909       direct or conformal conversion.  Direct conversion preserves radial
910       distance from the origin; conformal conversion preserves local angles,
911       so that each small-enough part of the image only appears to be scaled
912       and rotated, not stretched.  Conformal conversion puts the radius on a
913       logarithmic scale, so that scaling of the original image plane is
914       equivalent to a simple offset of the transformed image plane.
915
916       If you use three or more dimensions, the higher dimensions are ignored,
917       yielding a conversion from Cartesian to cylindrical coordinates, which
918       is why there are two aliases for the same transform.  If you use higher
919       dimensionality than 2, you must manually specify the origin or you will
920       get dimension mismatch errors when you apply the transform.
921
922       Theta runs clockwise instead of the more usual counterclockwise; that
923       is to preserve the mirror sense of small structures.
924
925       OPTIONS:
926
927       d, direct, Direct
928          Generate (theta,r) coordinates out (this is the default);
929          incompatible with Conformal.  Theta is in radians, and the radial
930          coordinate is in the units of distance in the input plane.
931
932       r0, c, conformal, Conformal
933          If defined, this floating-point value causes t_radial to generate
934          (theta, ln(r/r0)) coordinates out.  Theta is in radians, and the
935          radial coordinate varies by 1 for each e-folding of the r0-scaled
936          distance from the input origin.  The logarithmic scaling is useful
937          for viewing both large and small things at the same time, and for
938          keeping shapes of small things preserved in the image.
939
940       o, origin, Origin [default (0,0,0)]
941          This is the origin of the expansion.  Pass in a PDL or an array ref.
942
943       u, unit, Unit [default 'radians']
944          This is the angular unit to be used for the azimuth.
945
946       EXAMPLES
947
948       These examples do transformations back into the same size image as they
949       started from; by suitable use of the "transform" option to unmap you
950       can send them to any size array you like.
951
952       Examine radial structure in M51: Here, we scale the output to stretch
953       2*pi radians out to the full image width in the horizontal direction,
954       and to stretch 1 radius out to a diameter in the vertical direction.
955
956         $x = rfits('m51.fits');
957         $ts = t_linear(s => [250/2.0/3.14159, 2]); # Scale to fill orig. image
958         $tu = t_radial(o => [130,130]);            # Expand around galactic core
959         $y = $x->map($ts x $tu);
960
961       Examine radial structure in M51 (conformal): Here, we scale the output
962       to stretch 2*pi radians out to the full image width in the horizontal
963       direction, and scale the vertical direction by the exact same amount to
964       preserve conformality of the operation.  Notice that each piece of the
965       image looks "natural" -- only scaled and not stretched.
966
967         $x = rfits('m51.fits')
968         $ts = t_linear(s=> 250/2.0/3.14159);  # Note scalar (heh) scale.
969         $tu = t_radial(o=> [130,130], r0=>5); # 5 pix. radius -> bottom of image
970         $y = $ts->compose($tu)->unmap($x);
971
972   t_quadratic
973         $t = t_quadratic(<options>);
974
975       Quadratic scaling -- cylindrical pincushion (n-d; with inverse)
976
977       Quadratic scaling emulates pincushion in a cylindrical optical system:
978       separate quadratic scaling is applied to each axis.  You can apply
979       separate distortion along any of the principal axes.  If you want
980       different axes, use t_wrap and t_linear to rotate them to the correct
981       angle.  The scaling options may be scalars or vectors; if they are
982       scalars then the expansion is isotropic.
983
984       The formula for the expansion is:
985
986           f(a) = ( <a> + <strength> * a^2/<L_0> ) / (abs(<strength>) + 1)
987
988       where <strength> is a scaling coefficient and <L_0> is a fundamental
989       length scale.   Negative values of <strength> result in a pincushion
990       contraction.
991
992       Note that, because quadratic scaling does not have a strict inverse for
993       coordinate systems that cross the origin, we cheat slightly and use $x
994       * abs($x)  rather than $x**2.  This does the Right thing for pincushion
995       and barrel distortion, but means that t_quadratic does not behave
996       exactly like t_cubic with a null cubic strength coefficient.
997
998       OPTIONS
999
1000       o,origin,Origin
1001          The origin of the pincushion. (default is the, er, origin).
1002
1003       l,l0,length,Length,r0
1004          The fundamental scale of the transformation -- the radius that
1005          remains unchanged.  (default=1)
1006
1007       s,str,strength,Strength
1008          The relative strength of the pincushion. (default = 0.1)
1009
1010       honest (default=0)
1011          Sets whether this is a true quadratic coordinate transform.  The
1012          more common form is pincushion or cylindrical distortion, which
1013          switches branches of the square root at the origin (for symmetric
1014          expansion).  Setting honest to a non-false value forces true
1015          quadratic behavior, which is not mirror-symmetric about the origin.
1016
1017       d, dim, dims, Dims
1018          The number of dimensions to quadratically scale (default is the
1019          dimensionality of your input vectors)
1020
1021   t_cubic
1022         $t = t_cubic(<options>);
1023
1024       Cubic scaling - cubic pincushion (n-d; with inverse)
1025
1026       Cubic scaling is a generalization of t_quadratic to a purely cubic
1027       expansion.
1028
1029       The formula for the expansion is:
1030
1031           f(a) = ( a' + st * a'^3/L_0^2 ) / (1 + abs(st)) + origin
1032
1033       where a'=(a-origin).  That is a simple pincushion expansion/contraction
1034       that is fixed at a distance of L_0 from the origin.
1035
1036       Because there is no quadratic term the result is always invertible with
1037       one real root, and there is no mucking about with complex numbers or
1038       multivalued solutions.
1039
1040       OPTIONS
1041
1042       o,origin,Origin
1043          The origin of the pincushion. (default is the, er, origin).
1044
1045       l,l0,length,Length,r0
1046          The fundamental scale of the transformation -- the radius that
1047          remains unchanged.  (default=1)
1048
1049       d, dim, dims, Dims
1050          The number of dimensions to treat (default is the dimensionality of
1051          your input vectors)
1052
1053   t_quartic
1054         $t = t_quartic(<options>);
1055
1056       Quartic scaling -- cylindrical pincushion (n-d; with inverse)
1057
1058       Quartic scaling is a generalization of t_quadratic to a quartic
1059       expansion.  Only even powers of the input coordinates are retained, and
1060       (as with t_quadratic) sign is preserved, making it an odd function
1061       although a true quartic transformation would be an even function.
1062
1063       You can apply separate distortion along any of the principal axes.  If
1064       you want different axes, use t_wrap and t_linear to rotate them to the
1065       correct angle.  The scaling options may be scalars or vectors; if they
1066       are scalars then the expansion is isotropic.
1067
1068       The formula for the expansion is:
1069
1070           f(a) = ( <a> + <strength> * a^2/<L_0> ) / (abs(<strength>) + 1)
1071
1072       where <strength> is a scaling coefficient and <L_0> is a fundamental
1073       length scale.   Negative values of <strength> result in a pincushion
1074       contraction.
1075
1076       Note that, because quadratic scaling does not have a strict inverse for
1077       coordinate systems that cross the origin, we cheat slightly and use $x
1078       * abs($x)  rather than $x**2.  This does the Right thing for pincushion
1079       and barrel distortion, but means that t_quadratic does not behave
1080       exactly like t_cubic with a null cubic strength coefficient.
1081
1082       OPTIONS
1083
1084       o,origin,Origin
1085          The origin of the pincushion. (default is the, er, origin).
1086
1087       l,l0,length,Length,r0
1088          The fundamental scale of the transformation -- the radius that
1089          remains unchanged.  (default=1)
1090
1091       s,str,strength,Strength
1092          The relative strength of the pincushion. (default = 0.1)
1093
1094       honest (default=0)
1095          Sets whether this is a true quadratic coordinate transform.  The
1096          more common form is pincushion or cylindrical distortion, which
1097          switches branches of the square root at the origin (for symmetric
1098          expansion).  Setting honest to a non-false value forces true
1099          quadratic behavior, which is not mirror-symmetric about the origin.
1100
1101       d, dim, dims, Dims
1102          The number of dimensions to quadratically scale (default is the
1103          dimensionality of your input vectors)
1104
1105   t_spherical
1106           $t = t_spherical(<options>);
1107
1108       Convert Cartesian to spherical coordinates.  (3-D; with inverse)
1109
1110       Convert 3-D Cartesian to spherical (theta, phi, r) coordinates.  Theta
1111       is longitude, centered on 0, and phi is latitude, also centered on 0.
1112       Unless you specify Euler angles, the pole points in the +Z direction
1113       and the prime meridian is in the +X direction.  The default is for
1114       theta and phi to be in radians; you can select degrees if you want
1115       them.
1116
1117       Just as the t_radial 2-D transform acts like a 3-D cylindrical
1118       transform by ignoring third and higher dimensions, Spherical acts like
1119       a hypercylindrical transform in four (or higher) dimensions.  Also as
1120       with t_radial, you must manually specify the origin if you want to use
1121       more dimensions than 3.
1122
1123       To deal with latitude & longitude on the surface of a sphere (rather
1124       than full 3-D coordinates), see t_unit_sphere.
1125
1126       OPTIONS:
1127
1128       o, origin, Origin [default (0,0,0)]
1129          This is the Cartesian origin of the spherical expansion.  Pass in a
1130          PDL or an array ref.
1131
1132       e, euler, Euler [default (0,0,0)]
1133          This is a 3-vector containing Euler angles to change the angle of
1134          the pole and ordinate.  The first two numbers are the (theta, phi)
1135          angles of the pole in a (+Z,+X) spherical expansion, and the last is
1136          the angle that the new prime meridian makes with the meridian of a
1137          simply tilted sphere.  This is implemented by composing the output
1138          transform with a PDL::Transform::Linear object.
1139
1140       u, unit, Unit (default radians)
1141          This option sets the angular unit to be used.  Acceptable values are
1142          "degrees","radians", or reasonable substrings thereof (e.g. "deg",
1143          and "rad", but "d" and "r" are deprecated).  Once genuine unit
1144          processing comes online (a la Math::Units) any angular unit should
1145          be OK.
1146
1147   t_projective
1148           $t = t_projective(<options>);
1149
1150       Projective transformation
1151
1152       Projective transforms are simple quadratic, quasi-linear
1153       transformations.  They are the simplest transformation that can
1154       continuously warp an image plane so that four arbitrarily chosen points
1155       exactly map to four other arbitrarily chosen points.  They have the
1156       property that straight lines remain straight after transformation.
1157
1158       You can specify your projective transformation directly in homogeneous
1159       coordinates, or (in 2 dimensions only) as a set of four unique points
1160       that are mapped one to the other by the transformation.
1161
1162       Projective transforms are quasi-linear because they are most easily
1163       described as a linear transformation in homogeneous coordinates (e.g.
1164       (x',y',w) where w is a normalization factor: x = x'/w, etc.).  In those
1165       coordinates, an N-D projective transformation is represented as simple
1166       multiplication of an N+1-vector by an N+1 x N+1 matrix, whose lower-
1167       right corner value is 1.
1168
1169       If the bottom row of the matrix consists of all zeroes, then the
1170       transformation reduces to a linear affine transformation (as in
1171       t_linear).
1172
1173       If the bottom row of the matrix contains nonzero elements, then the
1174       transformed x,y,z,etc. coordinates are related to the original
1175       coordinates by a quadratic polynomial, because the normalization factor
1176       'w' allows a second factor of x,y, and/or z to enter the equations.
1177
1178       OPTIONS:
1179
1180       m, mat, matrix, Matrix
1181          If specified, this is the homogeneous-coordinate matrix to use.  It
1182          must be N+1 x N+1, for an N-dimensional transformation.
1183
1184       p, point, points, Points
1185          If specified, this is the set of four points that should be mapped
1186          one to the other.  The homogeneous-coordinate matrix is calculated
1187          from them.  You should feed in a 2x2x4 PDL, where the 0 dimension
1188          runs over coordinate, the 1 dimension runs between input and output,
1189          and the 2 dimension runs over point.  For example, specifying
1190
1191            p=> pdl([ [[0,1],[0,1]], [[5,9],[5,8]], [[9,4],[9,3]], [[0,0],[0,0]] ])
1192
1193          maps the origin and the point (0,1) to themselves, the point (5,9)
1194          to (5,8), and the point (9,4) to (9,3).
1195
1196          This is similar to the behavior of fitwarp2d with a quadratic
1197          polynomial.
1198

AUTHOR

1200       Copyright 2002, 2003 Craig DeForest.  There is no warranty.  You are
1201       allowed to redistribute this software under certain conditions.  For
1202       details, see the file COPYING in the PDL distribution.  If this file is
1203       separated from the PDL distribution, the copyright notice should be
1204       included in the file.
1205
1206
1207
1208perl v5.32.0                      2020-09-17                      Transform(3)
Impressum