
<html><body><div style="font-family: arial, helvetica, sans-serif; font-size: 12pt; color: #000000"><div>That works great! Thanks Evan!</div><div><br></div><hr id="zwchr" data-marker="__DIVIDER__"><div data-marker="__HEADERS__"><b>From: </b>"Even Rouault" <even.rouault@spatialys.com><br><b>To: </b>gdal-dev@lists.osgeo.org<br><b>Cc: </b>"Hays Barrett" <hbarrett@edac.unm.edu><br><b>Sent: </b>Tuesday, April 21, 2020 7:18:04 AM<br><b>Subject: </b>Re: [gdal-dev] Help with finding residual error for a given GCP.<br></div><div><br></div><div data-marker="__QUOTED_TEXT__"><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">Hays,</p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">GCPsToGeoTransform() is not used directly by gdalwarp (although it should likely be equivalent to running -order 1 polynomial), but by other code (OZI .map exporting, .tab exporting)</p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">You rather want to use the gdal.Transformer() API which is tested in</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">https://github.com/OSGeo/gdal/blob/master/autotest/gcore/transformer.py</p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">Even</p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> Hello all,</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> I want to get the error for each GCP before warping so I can calculate the</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> RMSE before georeferencing. With a list of GCPs, is it possible to get the</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> eventual map location of each pixel coordinate to calculate the</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> transformation error for each GCP? I have tried to use GCPsToGeoTransform</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> but that gave me a map location that was further from the intended chosen</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> GCP coordinate than the actual post- georeferencing pixel location.</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> Here is my GCPsToGeoTransform code:</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> ##CODE_START</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> from osgeo import gdal</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> #list of GCPs</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> gcpjson={'0': {'row': 5305.626973704321, 'col': 1158.432032601927, 'lat':</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> -107.18221364936802, 'lon': 35.08107313264506}, '1': {'row':</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> 778.3230540174817, 'col': 2317.0737231100225, 'lat': -107.17634116783263,</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> 'lon': 35.09958594602314}, '2': {'row': 3950.8777704175645, 'col':</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> 4880.465723118029, 'lat': -107.163460244414, 'lon': 35.08679518563491},</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> '3': {'row': 7609.882307742781, 'col': 4407.027704042081, 'lat':</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> -107.16570574683065, 'lon': 35.07164256912238}, '4': {'row':</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> 1822.0904204049455, 'col': 6486.900377637105, 'lat': -107.15529960157271,</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> 'lon': 35.096105471961536}, '5': {'row': 6291.396402113635, 'col':</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> 7216.1033738731185, 'lat': -107.15180240538085, 'lon': 35.077317337370516}}</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> #New list to put in gdal.GCP objects into</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> gcpList=[]</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> #Load GCP objects into gcpList</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> for gcp in gcpjson:</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> gcpList.append(gdal.GCP(gcpjson[gcp]['lat'],gcpjson[gcp]['lon'],0,gcpjson[gc</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> p]['col'],gcpjson[gcp]['row']))</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> #create GeoTranform from gcp list</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> GT = gdal.GCPsToGeoTransform(tuple(gcpList))</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> #Get pixel location from first GCP as a test.</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> col=gcpjson['0']['col']</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> row=gcpjson['0']['row']</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> xp = GT[0] + col*GT[1] + row*GT[2]</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> yp = GT[3] + col*GT[4] + row*GT[5]</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> print("POINT ("+str(xp)+" "+str(yp)+" 0)")</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> ##CODE_END</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> This prints out the WKT 'POINT (-107.18216298 35.0808978834 0)' which is</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> kinda close to -107.18221364936802, 35.08107313264506 but looking at my</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> post warped image the chosen pixel is right on top of the intended</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> coordinate so I am clearly doing something (or a lot) wrong.</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">></p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> Any help would be great.</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> Thanks,</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">> -Hays</p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;" data-mce-style="-qt-paragraph-type: empty; -qt-block-indent: 0; text-indent: 0px; margin: 0px;"> </p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">--</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">Spatialys - Geospatial professional services</p><p style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;" data-mce-style="-qt-block-indent: 0; text-indent: 0px; -qt-user-state: 0; margin: 0px;">http://www.spatialys.com</p><br></div></div></body></html>