Sometimes we have a map stored as a computer image but we don’t know how to convert between (E,N) map coordinates (e.g. British national grid) and (X,Y) image coordinates. For instance, we might require
The conversion can be found by choosing a set of easily-identified reference points, measuring their map coordinates and image coordinates, and applying the formulae below.
In theory the method will work if only 2 reference points are measured, but this doesn’t allow any check of the result. It’s best to measure several points, so that slight errors in measurement are averaged out, and gross errors can be found and corrected. The diagram shows an example with 6 reference points, labelled 0–5:
Suppose there are k reference points, with map coordinates
E[i], N[i] (i = 0, 1, ..., k-1)
and image coordinates
X[i], Y[i] (i = 0, 1, ..., k-1)
First calculate the mean coordinates, that is
Em = sum(E[i])/k; Nm = sum(N[i])/k; Xm = sum(X[i])/k; Ym = sum(Y[i])/k
Let lower-case letters denote coordinates relative to the mean values, that is
e[i] = E[i] - Em; n[i] = N[i] - Nm; x[i] = X[i] - Xm; y[i] = Y[i] - Ym
Calculate
P = sum(e[i]*x[i] - n[i]*y[i])
Q = sum(e[i]*y[i] + n[i]*x[i])
To find the conversion from map to image coordinates, calculate
R = sum(e[i]^2 + n[i]^2)
a = P/R
b = Q/R
c = Xm - a*Em - b*Nm
d = Ym - b*Em + a*Nm
Then the map-to-image conversion takes (E,N) to (X,Y), where
X = a*E + b*N + c
Y = b*E - a*N + d
The conversion from image to map coordinates is found in the same way, but with map and image coordinates exchanged. Calculate
R = sum(x[i]^2 + y[i]^2)
a = P/R
b = Q/R
c = Em - a*Xm - b*Ym
d = Nm - b*Xm + a*Ym
Then the image-to-map conversion takes (X,Y) to (E,N), where
E = a*X + b*Y + c
N = b*X - a*Y + d
It’s a good idea to check the conversion by applying it to the reference points. E.g. after working out the map-to-image conversion, apply it to each reference point (E[i],N[i]) and check that the result equals (X[i],Y[i]), within a reasonable tolerance to allow for errors in measuring. A large error means that the coordinate values used are incorrect. (If there are only 2 reference points this check is not available, because in this case the method gives a perfect fit for any input values.)
The algorithm does a least-squares fit over the set of transformations generated by rotation, scaling, and translation. Thus it assumes that shapes are the same in the image as in the map. In theory shapes will vary slightly according to the map projection used in the image, but in view of the low resolution of computer images the differences are likely to be negligible.