Landon Yarrington <lyarring...@gmail.com> writes: > First, the three control extension points listed in the map info are quite > far away and not on the map. And, the coordinates given for the points do > not appear to align with the actual locations when viewed on Google Earth. > Maybe this is a clue to the offset?
This is exactly a clue. The lat/lon given in the map data are labeled astronomic, and they are the results of some number of measurements which have considerable error. They decided to fix those stations with those values; this is the definition of a datum. In the 1940s I am aware of no better way to connect coordinates in Haiti to coordinates in CONUS. Surely the Army Map Service was world class in doing things like that. > Fort Islet lighthouse 18°33'31.33" -72°20'59.03", should be 18°33'18.08" > > -72°21'04.50" > Santo Domingo lighthouse 18°27'53.64" -69°52'52.74", should be 18°27'51.39" > -69°52'34.90" > Dame Marie 18°36'47" -74°25'53", unknown but possibly 18°36'57.21" > -74°25'26.94" You say "should be", but you are really saying "In some particular modern coordinate system". The definition of 0 longitude even in ITRF is totally arbitrary. It is just that the values adopted for mapping Haiti and the values measured in a modern frame are different. Yes, for latitude, one can say they are wrong, but gravity variation makes this hard. Looking only at minutes, and local - wgs84 lat lon Fort Islet 13" 5.5" Santo Domingo 2.2" -18" Dame Marie -10" -26" so this is not really that consistent. This leads me to think there is rotation too. But really, you have a map which you hope is internally consistent and you can identify points on the modern ground that match points on the map. > I took your advice and used the QGIS georeferencer. This is actually what I > initially did before reaching out here. I put google and bing satelite > imagery on the canvas to pin ground control points. I chose 19 GCPs and set There is the question of whether those imagery layers are aligned to recent WGS84, recent NAD83, something else, or are just off. But probably they are decent. I have checked MassGIS imagery against NAD83(2011) epoch 2010.0 (which the imagery is said to be in) and it matches to a pixel. But it was taken recently, in an area where there is great ground control and a statewide RTK network. > the algo to thin plate spline. The result has excellent accuracy. However, > the result is not easily reproducible. I have personally been to all the > locations I chose as GCPs and I know that---for the purpose of > georeferencing---those sites are unchanged in any meaningful way since the > 1942 aerial photography used to make the map was taken. It's great you have enough local knowledge to know what's undisturbed, and this is basically necessary to do what you are doing correctly. > Not everyone is going to have this specialized knowledge, and why > should they when there is already a UTM grid on the map? Can the grid > intersections make the georeferencing process more reproducible > without the need of ethnographic fieldwork? I think you are fundamentally misunderstanding the nature of UTM. UTM is a projection, not a datum. It takes latitude and longitude -- in some datum, relative to something -- and turns it into easting and northing. This is just math, and has no bearing on or input from whether the latitude/longitude pairs are consistent with modern WGS84. Further, this is NAD27-style UTM with the Clarke ellipsoid, not NAD83/WGS84-style UTM with the GRS80 ellipsoid. Assuming that the UTM coordinates are in NAD27 is exactly as wrong an assumption as assuming that the map lat/lon are in NAD27. > The map in question is the 1st edition ( > https://www.loc.gov/resource/g4940m.gct00206/?sp=168). There are at least 4 > subsequent editions (here is the 3rd > https://www.loc.gov/resource/g4940m.gct00206/?sp=167). The 1st edition has > the map info I quoted in my first message but is dropped from all later > editions. All later editions explicitly state the datum is NAD27 (or NAD83 > for 4th and 5th editions). The link you sent is the 2nd, but that's better, as it is the first that says NAD27, and the compilation date is 1963. I am not aware of any technology that would enable extending NAD27 from CONUS to Haiti in 1963. That seems early for TRANSIT, but the dates in wikipedia do not 100% preclude AMS trying to tie in Haiti with it: https://en.wikipedia.org/wiki/Transit_(satellite) If anyone knows any other techniques for tying Florida to Haiti in 1963, please post. I know this isn't really Surveying-L but I bet a lot of readers would be interested. > Curiously, the geographic representation in the > 1st edition is more expansive than the others; even though all the editions > have their neatlines within 19°50', -72°15' and 19°40', -72°00', the 1st > edition shows considerably more area. Example: the town of Limonade at the > bottom center of the map. The road that leads southeast from Limonade > extends about 1.29 km further than what is shown in later editions. In > fact, if you look at the grid intersection 2178000m N and 801000m E on the > 1st edition and the same point on the later editions, the grid is nearly > 1000 meters off. Compare > https://www.loc.gov/resource/g4940m.gct00206/?sp=168&r=0.314,0.573,0.372,0.229,0 > to > https://www.loc.gov/resource/g4940m.gct00206/?sp=167&r=0.284,0.495,0.238,0.147,0. > What is going on here?? Is this a datum difference, projection difference, > both?? If the neatline are the same lat/lon and they show different content, then either they have fixed errors or the datum is different. Same with projected coordinates. But the datum is different; they say so. I would get the NAD27 coordinates for their 3 control points and then compute distance/bearing from "AMS Haiti Datum of 1948" (I'm making that up) to NAD27. Those aren't real distances, but a datum shift. > So... the second thing I tried was selecting GCPs at each grid intersection > on the 1st edition but subtracting 1000 meters from the northings. (I did > 141 GCPs.) This got the map's accuracy much better, but there was still a > systematic offset. This doesn't really make sense. You are assuming that the original coordinate system is the same as some modern one and it just isn't. Think about if you were given a map in NAD27 and you are trying to georeference it into WGS84, but nobody will given you a transform between the two datums. Or how about you had an army of surveyors to do an experiment. Go to modern place that is 20km by 20 km. Pick 3 places and drive pins. For each, get really good GPS-derived coordinates. Then, generate 6 random numbers and add them into the 3 pairs of coordinates. Write those down and hand them to the surveying teams. The teams get to use total stations, with 1 arc second theodolites, and 1-2 ppm electronic distance measurement, but are forbidden from using any satellites, any outside info, and from making any astronomic observations. They are to survey the town. The map they make will be highly consistent with itself. But there will be some offset and rotation because the assumed coordinates of the 3 control points are in error. This is, as I understand it, essentially what happened, except that there was no oracle with random numbers. There was instead the errors in astronomic determination of latitude and longitude. Besides measuremnt errors, the direction of gravity is not aligned with the normal to the ellipsoid, and this results in latitude errors. > The best accuracy was in the center of the map. In the > 4th edition of the map, there is a note instructing how to convert from > NAD27 to NAD83, which reads, "Grid: Add 35mE, add 198mN. Geographic: > Subtract 1.4" Long, Add 2.2" Lat." OK. So, I also did these shifts to the > GCPs that I had already offset by 1000m N. For each easting, I added 35m > and each northing I added 198m. Then, I ran gdalwarp (code below). The > result still has an offset but it's the best I've been able to get using > the UTM grid, but sadly it is not nearly as accurate as the QGIS visual > approach. You are just fudging numbers by treating the UTM grid as true, and basically computing a single E/N shift by iteratively adjusting and looking. But you are not estimating a rotation. When you use GCPs with map coordinates in pixels and accurate WGS84 from imagery (if that is indeed acccurate), the georeferencer is doing least squares among all the points. > Again, I am indebted for any help. The code below only uses 8 GCPs of the > UTM intersections offset and shifted described above. > > ``` > # download the tiff from Library of Congress > wget > https://www.loc.gov/resource/g4940m.gct00206/?download=https%3A%2F%2Ftile.loc.gov%2Fstorage-services%2Fmaster%2Fgmd%2Fgmd4m%2Fg4940m%2Fg4940m%2Fgct00206%2Fcs000168.tif > -O cs000168.tif > # add 8 ground control points using the UTM grid described above > gdal_translate -of GTiff -gcp 833.408 5291.81 789035 2179198 -gcp 1146.904 > 5296.971 790035 2179198 -gcp 3963.057 6290.092 799035 2176198 -gcp 843.374 > 4660.743 789035 2181198 -gcp 1763.341 5935.218 792035 2177198 -gcp 2077.227 > 5941.007 793035 2177198 -gcp 8369.034 5735.588 813035 2178198 -gcp 8684.098 > 5741.097 814035 2178198 cs000168.tif cs000168_GCPs.tif > # warp the map with thin plate spline > gdalwarp -r near -tps -co COMPRESS=JPEG -co TILED=YES -s_srs EPSG:32618 > -t_srs EPSG:32618 cs000168_GCPs.tif cs000168_TO_epsg32618.tif > ``` Why are you running gdal_translate rather than just using the georeferencer and letting it make a new tif? Another consideration is that the map datum is a local one, but presumably distances are correct. So there should only be one translation and one rotation to go from map datum to NAD27. Or, one translation, one rotation, and one scale factor to go from map datum UTM to WGS84 UTM. And, the scale factor can be computed since the translation is small. By using a complicated mapping fit, you are allowing distortions in the map datum to be leveled out. If the map is distorted locally because of survey errors, that's good. But if the map is internally consistent and the individual points have random errors, that's bad. _______________________________________________ PROJ mailing list PROJ@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/proj