5

I am trying to georeference an image using many (>1,000,000) GCPs (derived by matching two images using triangulation). However, after many attempts doing this in different ways I think I found that there is a limit to the number of GCPs that I can use.

ds = gdal.Open("tmp_GCP.TIF", gdal.GA_Update)
wkt = ds.GetProjection()
ds.SetGCPs(gcp_list_gdal[:10923], wkt)
print(ds.GetGCPCount()) # output 10923
ds.FlushCache()
ds = None
gdal.Open("tmp_GCP.TIF").GetGCPCount() # output 0

When I set 10923 GCPs I obtain a file that seemingly has zero GCPs and when I set this to 10922 I do get a file with 10922 GCPs. This seems to be independent of the width/height of the raster.

Also, I tried to set the GDAL cache to a higher value, but this also didn't have an effect.

gdal.SetConfigOption('GDAL_CACHEMAX', '1000000')

Edit: Just found out that I am getting the following error in the terminal when running this code in a Jupyter Notebook (which I don't get when I use gcp_list_gdal[:10922]).

Warning 1: TIFFFetchNormalTag:Incorrect count for "GeoTiePoints"; tag ignored
Warning 1: TIFFFetchNormalTag:Incorrect count for "GeoTiePoints"; tag ignored

I have tried this on Ubuntu 16 and Ubuntu 18, using Gdal 2.2, 2.4 and 3.0, but I always get the same error.

Edit 2: Because of the error I thought it might be a problem with the amount of GCP information that can be stored in a GeoTIFF file, so I tried another approach in python (trying not to store the DS before applying a warping), but it seems the same problem occurs (since it does produce an output using 10922 GCPs but not when using 10923 GCPs).

ds = gdal.Open('tmp_GCP.TIF', gdal.GA_Update) 
gcp_srs = osr.SpatialReference()
gcp_srs.ImportFromEPSG(4326)
wkt = gcp_srs.ExportToWkt()
ds.SetGCPs(gcp_list_gdal[:10922], wkt)
dstDS = gdal.Warp('TEST.TIF', ds, format='GTiff', tps=True, dstNodata=0,
                  srcNodata=0, resampleAlg=gdal.GRA_Bilinear)

Edit 3: I found out about the option to use geolocation arrays https://gdal.org/development/rfc/rfc4_geolocate.html which seems very interesting. As I understand it you can basically create two bands of the same size as the image and indicate for each pixel the desired longitude and latitude of the pixel. This should work better with a large number of GCPs so I am going to try this.

3
  • 1
    Out of curiosity, why do you need >1,000,000 GCPs? Nov 29, 2019 at 11:18
  • We want our image to match as much as possible with our calculated height map (calculated using stereo processing).
    – David
    Nov 29, 2019 at 11:46
  • This is actually related to a previous question I posted about not being able to orthorectify my WorldView-3 images using the RPC model (gis.stackexchange.com/questions/328366/…). So that's why I am trying to use GCP points directly.
    – David
    Nov 29, 2019 at 12:08

1 Answer 1

1

The error seems to come from the libtiff library, and the corresponding source code file is https://gitlab.com/libtiff/libtiff/-/blob/master/libtiff/tif_dirread.c

So the answer to the question "Is there a limit to the number of GCPs that can be used?" is yes, there is a limit. However, the limit is maybe not set to number of GCPs but to the size of tags in bytes https://gitlab.com/libtiff/libtiff/-/blob/master/libtiff/tif_dirread.c#L1260

/* Caution: if raising that value, make sure int32 / uint32 overflows can't
 * occur elsewhere */
#define MAX_SIZE_TAG_DATA 2147483647U

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.