Dealing with Coordinates: Coordinate Transformation From Source to Map and back

1
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
YS
YM
P
P
XM(P)
XS(P)
YM(P)
YS(P)
XM
XS
2
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Map
Source
P1
YS(P1)
DS
P2
dYS
YS(P2)
dXS
XS(P2)
XS(P1)
dXS XS(P2) - XS(P1) dYS YS(P2) - YS(P1) DS
Sqrt(dXS2 dYS2)
dXM XM(P2) - XM(P1) dYM YM(P2) - YM(P1) DM
Sqrt(dXM2 dYM2)
Scaling Scale DM / DS dXM / dXS dYM /
dYS see note below Source to Map XM(P)
Scale XS(P) YM(P) Scale YS(P) Map to
Source XS(P) XM(P) / Scale YS(P)
YM(P) / Scale
Note The terms in ... are only valid if
no Rotations applies
3
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Coordinate Shift
YM
Assume Source coordinates have the same scale
as Map coordinates and Source is overlayed over
Map such that corresponding points (eg. P1, P2,
etc) match up. Visualise Source as scaled aerial
photo on transparent film fitted to
base-map. Because reference is made to different
origins in both, we can observe that source
origin appears to coincide with a map-position
different from map-origin (or visa versa). In
order to translate source coordinates for
points into those of the map-coordinate system,
we must apply a shift in X and Y direction each.
The amount of shift is the negative of the value
of one systems origin in coordinates of the
other. Using the map origin, its position in the
source system is PMO (XS(MO), YS(MO))
YS(P1)
YS(P2)
XM
YS(MO)
XS(MO) XS(P1) XS(P2)
Source to Map XM(P) XS(P) - XS(MO) and
YM(P) YS(P) - YS(MO) Map to Source XS(P)
XM(P) XS(MO) and YS(P) YM(P)
YS(MO) Following the above instruction, the
Source Origin expressed in map coordinates
is PSO (XM(SO), YM(SO)) (-XS(MO),
-YS(MO)) Map to Source XS(P) XM(P) - XM(SO)
and YS(P) YM(P) - YM(SO)
Same Structure
4
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Rotating Coordinate Systems
Assume Source coordinates have the same scale
as Map coordinates and Source is overlayed over
Map such that corresponding points (eg. P1, P2,
etc) match up. Visualise Source as scaled aerial
photo on transparent film fitted to base-map.
Assume further, the origins of both systems also
coincide, but the angular orientation of the
coordinate axis differs by angle a (alpha).
YM
YS
a
P1
XS
P2
Thus, the coordinates for identical points (eg.
P1) differ from system to system. In order to
transform coordinates from the source system into
those of the map system we need to formulate the
following relationship XM(P) xa - xb
where xa XS(P) cos(a) and
xb YS(P) sin(a) XM(P) XS(P) cos(a) -
YS(P) sin(a) and similarly YM(P) XS(P) sin(a)
YS(P) cos(a)
YS(P1)
yb
a
xb
XS(P1)
ya
XM
xa
XM(P1)
5
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Putting it all together So far we have assumed
that before we apply a new transformation type
(1) to (3), the previous step has already been
accomplished. For example, the XS, YS coordinate
values where assumed to have been already shifted
to make the origins of both systems coincide and
that before step (2) occurred, the source
coordinates were already scaled to match that of
the map. Thus we can combine these steps in
reverse order (3)-gt(2)-gt(1) into one unique set
of transformation equations accomplishing all
three steps at the same time Substitute Shift
into Rotation XM(P) (XS(P) - (XS(MO)) cos(a) -
(YS(P)-YS(MO) sin(a) YM(P) (XS(P) - (XS(MO))
sin(a) (YS(P)-YS(MO) cos(a)
Scaling XM(P) Scale XS(P) YM(P) Scale
YS(P)
(1)
Shift XM(P) XS(P) - XS(MO) YM(P) YS(P) -
YS(MO)
(2)
Rotation XM(P) XS(P) cos(a) - YS(P)
sin(a) YM(P) XS(P) sin(a) YS(P) cos(a)
(3)
Now we can apply the scaling to the right hand
sides as in (1) above and obtain the final
COORDINATE TRANSFORMATION EQUATIONS XM(P)
Scale (XS(P) - (XS(MO)) cos(a) - Scale
(YS(P)-YS(MO) sin(a) YM(P) Scale (XS(P) -
(XS(MO)) sin(a) Scale (YS(P)-YS(MO) cos(a)
6
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Exploring the COORDINATE TRANSFORMATION
EQUATIONS in shorthand
XM(P) Scale (XS(P) - (XS(MO)) cos(a) - Scale
(YS(P)-YS(MO) sin(a) XNew X0 a
XOld - b YOld YM(P) Scale (XS(P) - (XS(MO))
sin(a) Scale (YS(P)-YS(MO) cos(a) YNew
Y0 b XOld a Yold
Where X0 -ScaleXS(MO)cos(a)
ScaleYS(MO)sin(a) Y0 -ScaleXS(MO)sin(a)
-ScaleYS(MO)cos(a) a Scale cos(a)
and b Scale sin(a) even shorter
shorthand X X0 a X - b Y
Y Y0 b X a Y for a Similarity Transform
Coordinates of points given in the Source
System are (1) shifted to a new Origin matching
the Map-Origin, (2) scaled to match the
Map-Scale and (3) subjected to a common
rotation a to match the orientation of the
content of the Map Thus, neither of these steps
will alter the geometric position
relationship between these points, the shape is
maintained, the resulting point constellation is
SIMILAR to what it was in the Source. Therefore
this kind of coordinate transformation is
called Similarity Transformation
7
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Some additional Thoughts on Coordinate
Transformations using the Similarity Transform as
Example
XM(P) Scale (XS(P) - (XS(MO)) cos(a) - Scale
(YS(P)-YS(MO) sin(a) X X0 a X - b
Y YM(P) Scale (XS(P) - (XS(MO)) sin(a)
Scale (YS(P)-YS(MO) cos(a) Y Y0 b
X a Y
In order to apply such a coordinate
transformation we need to know the amount of (1)
Shift by XS(MO) and YS(MO), (2) Scale DM / DS
(see slide 2) (3) Rotation angle a Though this
might be the case under special circumstance, it
is more likely that we are not supplied with this
data. In most cases (as in the example of using
an areal photo as Source and a cadastral map as
Base Map) we are required to determine these
transformation parameters from coordinates of
corresponding points in both systems. Recognising
also, that the four parameters are implicitly
contained in the short-form parameters X0, Y0, a
and b, and that for every point for which we know
coordinates in both systems a pair of X .....
and Y .... Equations can be written, we may
solve for the short form parameters by using two
such points resulting in four equations to
solve for four unknown parameters. In this
respect it is better to re-write the short-form
equations as
X(i) X0 a X(i) - b Y(i) Y(i) Y0 b X(i)
a Y(i)
Similarity Transform
where (i) is an index denoting the ith point
known in both systems. Once X0, Y0 and a and b
are computed the transformation can be applied to
all other points given in Source coordinates.
8
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Another look at Transformation Equations
XM(P) Scale (XS(P) - (XS(MO)) cos(a) - Scale
(YS(P)-YS(MO) sin(a) X(i) X0 a
X(i) - b Y(i) YM(P) Scale (XS(P) - (XS(MO))
sin(a) Scale (YS(P)-YS(MO) cos(a)
Y(i) Y0 b X(i) a Y(i)
The above form of the coordinate transformation
equations assumes a unique scale exists between
Source and Map. This is not always true. Again
take aerial photography as Source as an example.
A true overall scale only exists here if and only
if the photograph is truly vertical and the
terrain is completely flat and horizontal. This
of course rarely happens. Assume completely flat
and horizontal terrain, but the photography not
being exactly vertical. In this case the scale of
the Source (photo) diminishes towards the far
horizon, and the use of exact projective
geometry would be appropriate. However, we may
approximate this situation by applying two
different scales to the XS and YS coordinates of
the Source. With this assumption the
transformation equations need to be amended
XM(P) XScale (XS(P) - (XS(MO)) cos(a) -
YScale (YS(P)-YS(MO) sin(b) X(i)
X0 a X(i) - b Y(i) YM(P) XScale (XS(P) -
(XS(MO)) sin(a) YScale (YS(P)-YS(MO) cos(b)
Y(i) Y0 c X(i) d Y(i)
This also requires a re-write of the shorthand
form of the equations to contain also six (6)
transformation parameters compared to four in the
initial similarity form. The same effect will be
encountered when we allow different rotation
angles in respect to the X and Y axis of either
system. This also is a typical trait of oblique
photography, where the map coordinate grid would
appear as non- rectangular if it would show.
a Xscale cos(a) b Yscale sin(b) c
Xscale sin(a) d Yscale cos(b) X0 and Y0
also will change, but retain a bias-character.
9
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
This form of coordinate transformation is called
Affine Transformation
XM(P) XScale (XS(P) - (XS(MO)) cos(a) -
YScale (YS(P)-YS(MO) sin(b) X(i)
X0 a X(i) - b Y(i) YM(P) XScale (XS(P) -
(XS(MO)) sin(a) YScale (YS(P)-YS(MO) cos(b)
Y(i) Y0 c X(i) d Y(i)
The longhand form of the equations is called
parametric, because the transformation
parameters contained have physical meaning, such
as shift, scale, rotation and skew (a -
b) d Skewness d can be visualised as the
departure of either or both of the coordinate
systems from being rectangular
d
The shorthand form is called non-parametric
since its parameters are not connected directly
to any physical meaning. In order to
obtain physical parameters we have to perform
additional computations tan (a) c / a
tan(b) b / d a atan(c/a)
b atan(b/d) Xscale a / cos(a)
c / sin(a) Yscale b / sin(b) d /
cos(b) similarly X0 and Y0 can be
retrieved. In most cases the non-parametric
approach is sufficient.
P1
YM(P1)
P2
YM(P2)
XM(P1)
XM(P2)
Example of a skewed coordinate system.
10
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Looking at Coordinate Transformations in a
generalised way
X(i) X0 a X(i) - b Y(i) Y(i) Y0 c X(i)
d Y(i)
Affine Transform Six transformation Parameters
required, thus at least three Control Points
with coordinates known in both systems are
required to determine these parameters. Similarit
y Transform a d and b c gt four parameters
required two Control Points Often more than the
minimum number of Control Points are available
for use. In such a situation it is important (a)
to use at least those points which cover the full
extend of the Source Point area from extreme to
extreme and (b) remaining required points equally
spaced inside the extremes, but not forming a
straight line with any other two points
and (c) preferably to use all available points
(if the transformation procedure (program)
permits). Note using more than the minimum
required number of Control Points is
theoretically permissible, but in practise local
systematic and random (error) variations in scale
and position will lead to contradictions. To
resolve these local contradictions the
transformation needs to incorporate a Least
Squares Adjustment approach. Similarity
Transforms with Least Squares Adjustment (LS) are
known in some literature as Helmert Transform.
Affine Transforms usually incorporate LS.
Furthermore, the LS approach allows also the
computation of a measure of reliability in form
of a Standard Deviation of Fit between both
coordinate systems. The smaller the standard
deviation the better the agreement between Source
and Map.
11
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Looking at Coordinate Transformations even more
generalised
X(i) X0 a X(i) - b Y(i) Y(i) Y0 c X(i)
d Y(i)
In mathematical sense, the two above equations
can be characterised as polynomial equations
of 1st order, since they include a constant (X0
or Y0) and two linear relationship terms (eg. a
X(i)) each. In this respect there is no limit to
expand these polynomial equations into higher
order, such as quadratic , cubic or even higher.
A quadratic form for example is
X(i) X0 a X(i) b Y(i) e X(i)Y(i) f
X(i)X(i) g Y(i)Y(i) Y(i) Y0 c X(i) d
Y(i) h X(i)Y(i) i(X(i)X(i) j Y(i)Y(i)
This non-parametric set of equations can be
used to also model the local systematic
distortion in a spatial parabolic/hyperbolic
form. (It is a 2nd order polynomial). Higher
order polynomial equations can deal with even
more complex local distortions. However, the
draw-back is the increase in the number
of transformation parameters (in the above set to
12, requiring at least 6 Control Points for their
determination. In practice this kind of approach
allows the correction of systematic distortions
in the fit between Source and Map. Remember Maps
are the result of projecting a curved surface
(earth) onto a flat sheet of paper causing
systematic distortions! One could find a likeness
for both in visualising a rubber balloon as a
possible Source (or Map) to be flattened out.
Transformations of this kind, using higher order
mathematical equations are often included in what
is known as RUBBER-SHEETING.
12
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Looking at Coordinate Transformations even more
generalised
X(i) X0 a X(i) b Y(i) Y(i) Y0 c X(i)
d Y(i)
Note Parameter b has been changed from -b to
b !!!
For Similarity Transforms a d and
b -c including X0 and Y0, 4 unknown
parameters need to be determined. This
required at least 2 Control Points ( 2
equations for two Control Points 4 equations
for 4 unknowns ) For Affine Transforms a, b, c
and d will in most cases be all different from
each other including X0 and Y0, 6 unknown
parameters need to be determined. This
required at least 3 Control Points ( 2
equations for three Control Points 6 equations
for 6 unknowns )
Expanding into 3-Dimensions The structure of the
two equations lends itself easily for an
expansion into 3D
In the general 3D Affine Transformation case, 12
parameters need to be determined. Since we can
write three equations per point, (12 / 3) 4
Control Points are required.
X(i) X0 a X(i) b Y(i) e Z(i) Y(i)
Y0 c X(i) d Y(i) f Z(i) Z(i) Z0 g
X(i) h Y(i) k Z(i)
13
Dealing with Coordinates Coordinate
TransformationFrom Source to Map and back
Looking at Coordinate Transformations even more
generalised Expanding into Matrix Algebra
Algebraic Form Matrix Form
X X0 R x X 3D Form X(i)
X0 where X Y(i) X0 Y0
Z(i) Z0
a b e X(i) R c d
f X Y(i) g h k
Z(i)
X(i) X0 a X(i) b Y(i) e Z(i) Y(i)
Y0 c X(i) d Y(i) f Z(i) Z(i) Z0 g
X(i) h Y(i) k Z(i)
Identical !
2D Form X X0 R x X simply leave the
Z-element out for all Vectors and delete the
third column and row from Matrix R (found on
the right).
X is the result Vector of transformed
Coordinates X is the Vector of given Source
Coordinates X0 is the Vector of Map Origin
expressed in Source Coord. R is the scaled
Rotation Matrix
N-Dimensional Form Simply add coordinate
identifiers to each Vector to make up N
coordinate-axis and add equally as many rows and
columns to Matrix R
For every point X, given in the Source System, a
set of equations as shown above can be written
and solved to find the coresponding Target (Map)
coordinate set X .
14
Dealing with Coordinates Coordinate
TransformationUTM Zone to Zone Transformation
In AMG/MGA System
On the Ellipsoid
y1(P)
x1(P)
x2(P)
y2(P)
F(P)
P
F(P)
P
Equator
CM1
CM2
l(P)
dl 3 degree for AMG 6 degree wide Zones
dl 3 degree for AMG 6 degree wide Zones
Remember At the eastern and western edges of
AMG/MGA-Zones, an additional 1/2 degree overlap
is provided. Points in this overlap should be
provided with coordinates in both
Zone- Systems. The normal procedure is to
determine point coordinates in one of the two
Zones by measurement (eg. GPS) and to
transform these into the other overlapping Zone.
15
Dealing with Coordinates Coordinate
TransformationUTM Zone to Zone Transformation
Assume points are given by coordinates in one
AMG/MGA Zone System and are also required to be
known by coordinates in a neighbouring
Zone. Zone to Zone Transformation of UTM
coordinates (AMG/MGA for Australia) can be
accomplished in two ways (1) Transform given
UTM coordinates into Longitude and Latitude
(geodetic coordinates) and subsequently
transform these back into UTM coordinates (as
Easting and Northing) of the required
adjacent Zone. This is a sequence of
Re-Projection and Projection. (See following
pages) (2) Using direct Zone-to-Zone
Transformation equations. These kind of
transformation equations are in most cases
series developments in complex polynomial
form N iE a0 ib0 (a1 ib1)(y ix)
(a2 ib2)(y ix)(y ix) ..... The
real and imaginary part of these polynomials are
finally split apart to find the solution for
E (Easting) and N (Northing). Note Method (1)
can be used to transform map-coordinates of any
projection type into any other projection type,
as long as appropriate projection equations are
used. Method (2) usually is restricted to
transformations of coordinates between two
identical projections types, only differing in
the positioning of the projection surface (eg.
Two zones of AMG). Some transformation equations
of type (2) also allow for a transformation
between different datum definitions, but the same
projection type (eg. AGD66 to AGD84 or
GDA94). For both approaches, datum parameters
are required and must be known (eg. ADG84
parameters).