Georegistration: A Worked Example

We're going to work through making a quad sheet and georegistering it. The quad is going to be [90W 29N 89W 30N]. If you're wondering why that quad out of all of the possible quads, it's because I know that 90W 30N hits one of my favorite spots: New Orleans. Let's pick a scale for our quad to be 1:100000, and let's give ourselves one inch margins. We're going to use UTM zone 16 WGS 84 as our coordinate reference system, and I'm going to use GEOTRANS to do the calculations. The the table of longitudes (λ), latitudes (ϕ), eastings and northings is:

  λ  ϕ  E(m)    N(m)
-90 29 207729 3211697
-90 30 210590 3322576
-89 30 307085 3320469
-89 29 305179 3209635

The bounding eastings and northings are:

[207729 3209635 307085 3322576]

Taking the difference between the upper right and lower left values and scaling the result leaves us with a rectangle 0.99356 m by 1.12941 m. To put us into PDF world, we're going to convert the meters to points and add 144 points to the result to give us our one inch margins. The resulting media box is (rounded values):

[0 0 2960 3345]

The neatline falls an inch (72 points) inside the page media box at:
[72 72 2888 3273]
What we've done here is aligned the axes of the projected and page coordinate systems. Moreover, we constrain the scales to be constant regardless of direction. This let's us reduce the problem from solving:

|a c H||h|   |x|
|b d V||v| = |y|
|0 0 1||1|   |1|

to

|a 0 H||h|   |x|
|0 a V||v| = |y|
|0 0 1||1|   |1|

The solution is given by:

a = (x1 - x0)/(h1 - h0)
H = x0 - a h0
V = y0 - a v0

Using the values:

x0 = 207729
y0 = 3209635
x1 = 307085
h0 = v0 = 72
h1 = 2888

yields:

a = 35.28267
H = 205188.64
V = 3207094.8

In PDF notation, the full CTM would be:

[35.28267 0 0 35.28267 205188.64 3207094.8]

This matrix defines the relationship between a point on the page and a point in the projected coordinate system. It's inverse takes the projected coordinate system into the page coordinate system:

[0.028343 0.0 0.0 0.028343 -5815.5645 -90897.17]

Using this, we can discover where the corners of the map neatline will fall:

  λ  ϕ  E(m)    N(m)  H(pt) V(pt)  
-90 29 207729 3211697   72   132
-90 30 210590 3322576  153  3274
-89 30 307085 3320469 2888  3215
-89 29 305179 3209635 2834    74

In a later post, we'll draw this map...