Rectification Question

[PhilipGriffith]

I'm currently trying to automate the rectification process of georeferenced aerial imagery stored in SID format.
To that end I'm using Python to pull the Source and Target GCPs from an image's associated aux.xml file and rectifying the image using the Warp tool. So far things have gone swimmingly. However, I've discovered some images that were georeferenced about a year ago, and while they have a defined spatial reference, they don't seem to have any associated GCPs. In other words, when I look at their aux.xml files, I don't see any lines like this:

- <SourceGCPs xsi:type="typens:ArrayOfDouble">
  <Double>2.968341898744812</Double>
  <Double>6.2462671352726753</Double>

Despite this fact, when I load them into ArcMap they appear correctly georeferenced. In fact, I can rectify them using the Georeferencing Toolbar without a problem. This would be a viable solution, except that there are over 400 images requiring rectification. This isn't an automation question, but rather a more basic one:

How is ArcMap correctly georeferencing the images if they don't seem to have any GCPs?

The images have an associated world file (.sdwx) such as the one below:

0.83235784953486269
0.026282018489577653
0.025363891070187279
-0.83526446120022968
640231.33534762322
3429918.1939442093

They also contain data within the GeodataXform tag of the aux.xml with coefficient information such as this:

- <CoeffX xsi:type="typens:ArrayOfDouble">
  <Double>719580.22190854128</Double>
  <Double>495.2226517556241</Double>
  <Double>75.839059207416327</Double>
  </CoeffX>

But if I'm understanding things correctly, the coefficients alone aren't enough to rectify an image.
The images in question also have an associated sid.xml file that contains information like the following:

<spatRepInfo>
- <Georect>
- <cellGeo>
  <CellGeoCd Sync="TRUE" value="002" />
  </cellGeo>
  <numDims Sync="TRUE">2</numDims>
  <tranParaAv Sync="TRUE">1</tranParaAv>
  <chkPtAv Sync="TRUE">0</chkPtAv>
- <cornerPts>
  <pos Sync="TRUE">719579.880946 3360574.066594</pos>
  </cornerPts>
- <cornerPts>
  <pos Sync="TRUE">719579.880946 3365032.095767</pos>
  </cornerPts>
- <cornerPts>
  <pos Sync="TRUE">724849.400304 3365032.095767</pos>
  </cornerPts>
- <cornerPts>
  <pos Sync="TRUE">724849.400304 3360574.066594</pos>
  </cornerPts>
- <centerPt>
  <pos Sync="TRUE">722214.640625 3362803.081180</pos>
  </centerPt>

This is the closest thing I can find to anything like a GCP. Is this what the Georeferencing toolbar is using when it rectifies the image?
If so, how would you alter the data to use as an input for the Warp tool?

Thanks for any and all help. I'm mystified!
-Philip

[PhilipGriffith · Accepted solution]

After doing some thinking I finally realized my error: the coordinates specified in the .sid.xml file under the <Georect> tag describe the four points and the center of the image in relation to a specific projection, which is how ArcMap knew where to correctly display them.

Given this fact, the same coordinates can serve as both the Source and the Target GCPs when using the Warp tool.
While I don't know which points were originally chosen to establish where the corners of an image should be georeferenced, I can use the spatial information associated with the corners to correctly rectify the image.

This was only an issue with files that did not include a .sid.aux.xml file. I'm not sure which tool (or program version) was used to georeference the images, but with the associated .sid.xml file their rectification was still able to be automated.

All the best,
Philip

[MichaelLindgren]

Hi Phillip,

I am contacting you to inquire about your way forward working with rectifying TIFF images.  I am having some difficulty working with the .tif.aux.xml file as well and your question and self-answer are the closest things I have found to what I am looking for.  I would like to write a small function in R to deal with a standard TIFF file that has been georeferenced in ArcMap 10, but not rectified to create a new GeoTIFF, but I also cannot for the life of me figure out what the heck those GCP values are in the .tif.aux.xml file that is created when the georeferencing took place.  I am thinking that using a combination of gdal_translate / gdalwarp should produce the GeoTIFF in the right locations given the GCP's.  Question is, how do I create those GCP's from the coefficients that were given in that xml file. 

Ultimately, would you be interested in discussing further how I may be able to exploit the info in that XML and warp the images either in ArcMap or in GDAL / R. 

Many thanks in advance,

Michael

[PhilipGriffith]

Hi Michael,

Sorry for the lagged response. Hopefully you've figured this out by now! If not, I'm curious to know if your xml files have the <SourceGCP> and <TargetGCP> tags. If you have those tags, the first two coordinate values in the source tag are the X and Y values of the first Source GCP. The second two are the second XY coordinates, etc., like this:

X
Y
X
Y

It's exactly the same for the Target GCPs. If you can extract these and append them to a list, it won't take much to plug them into whatever warp tool you're using.

I didn't know you can warp images using R! I can run basic statistical models, but I'm still learning how to use its spatial toolsets. Anyway, I hope this gets you a little closer to a solution.

All the best,
Philip

[MichaelLindgren]

Thanks soo much for your response!  The xml files did have the needed tags and using the tools in the GDAL Python bindings I was able to use that information to rectify the images.  It basically grabs those tags from the xml, creates GCPs from them and rectifies the newly generated file using the gdalwarp command ( had to use this as there are not yet Python bindings to the warp api in GDAL ). 

Another way to solve this is to use gdal_translate to add the gcps to the image and then run gdalwarp on that new image with the GCPs to rectify it.

Further, using the absolutely amazing Raster package in R there is a command named 'rectify' that may do the trick if the gcps exist in the geotiff header, but that is an untested solution.

The code for the tool I wrote (tested only on the data needed for my particular end), is available for download at: https://github.com/EarthScientist/GDAL_tools/blob/master/gcpRectify.py

More documentation on the code to follow in the coming days as I have some time to work on it.

Hope someone finds it useful!

Michael Lindgren

[PhilipGriffith]

Glad I could help!