[gdal-dev] Possible concurrency issue in OGRCoordinateTransformation

[Wood, Alexander]

Hi,
I've encountered what appears to be a concurrency issue when calling OGRCoordinateTransformation::Transform from multiple threads.  What I've found is that when Transform is called from multiple asynchronous tasks, I'm getting inconsistent results given the same inputs.  This is pronounced with the projection used in the GINA Elevation dataset (Alaska); sample data can be found here:
http://ifsar.gina.alaska.edu/data/2013/DTM/IFSAR.SDMI.2013.FUGRO.DTM_N59W144/IFSAR.SDMI.2013.FUGRO.DTM_N59W144.tar.gz

The projection used by this datasource is:
PROJCS["NAD_1983_CORS96_Alaska_Albers",GEOGCS["GCS_NAD_1983_CORS96",DATUM["NAD_1983_CORS96",SPHEROID["GRS_1980",6378137,298.257222101]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],PROJECTION["Albers_Conic_Equal_Area"],PARAMETER["standard_parallel_1",55],PARAMETER["standard_parallel_2",65],PARAMETER["latitude_of_center",50],PARAMETER["longitude_of_center",-154],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]]]

Or in proj4:
+proj=aea +lat_1=55 +lat_2=65 +lat_0=50 +lon+0-154 +x_0=0 +y_0=0 +ellps=GRS80 + units=m +no_defs

The following code snippet reproduced the problem for me.  If I were comment out the async tasks and simply perform the transform synchronously, the outputs all remain consistent with the expected result.  Also make note of the CPLConfigOption "USE_PROJ_480_FEATURES".  Setting this to "NO" will also produce consistent results, as indicated in the inline comments below.

bool TestGdalConcurrency(const std::string& filename)
{
    auto dataset = GDALOpen(filename.c_str(), GDALAccess::GA_ReadOnly);
    auto projectionName = GDALGetProjectionRef(dataset);

    // By default, this is set to yes
    // Configuring this setting to "NO" will force GDAL to acquire a lock when calling into proj
    // This is a workaround, as we'd prefer lockfree concurrent access to proj.
    CPLSetConfigOption("USE_PROJ_480_FEATURES", "YES");

    // Setup "From" SpatialReference
    OGRSpatialReference wgs84SpatialReference;
    wgs84SpatialReference.SetWellKnownGeogCS("WGS84");
    // Setup "To" SpatialReference
    OGRSpatialReference spatialReference(projectionName);
    // Create transformation from WGS84 to the sourcefile's projection
    auto pFromPlanet = OGRCreateCoordinateTransformation(&wgs84SpatialReference, &spatialReference);

    // Will also provide this consistent latitude, longitude as input
    // The expectation is that given constant inputs, we should always receive the exact same output.
    double longitude = -143.8996;
    double latitude = 60.0074;
    double x = longitude;
    double y = latitude;

    // Our truth data (we're confident in this answer since there's no concurrent access into proj yet.
    pFromPlanet->Transform(1, &x, &y, nullptr);

    // Asynchronously call Transform from as many threads as I can get, every result should match the truth data
    std::vector<std::future<bool>> results;
    for (int i = 0; i < 100000; i++)
    {
        results.emplace_back(std::async(std::launch::async, [longitude, latitude, pFromPlanet, x, y]()
        {
            double _x = longitude;
            double _y = latitude;

            pFromPlanet->Transform(1, &_x, &_y, nullptr);
            return _x == x && _y == y;
        }));
    }

    // Sum the number of failed comparisons
    auto numFailures = 0;
    for (auto& result : results)
    {
        if (!result.get())
            ++numFailures;
    }

    OGRCoordinateTransformation::DestroyCT(pFromPlanet);
    GDALClose(dataset);

    return numFailures == 0;
}

If there's something I'm missing that would allow for lock-free coordinate transformations I'd be happy to hear about it.

Thanks,
Alex
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20160913/6f3e11a5/attachment-0001.html>