Skip to content
Permalink
Browse files

Update qgsrasterfilewriter.cpp

Bad geo-referencing bug fix

Using VRT with large coordinate system (like RGF93/CCxx) and high
resolution raster, highest resolution tiles were shifted regardingi
original raster.
It seemed that when creating the VRT file, georeferencing is
truncated while converted to text, leading to loss of precision.
  • Loading branch information
gui2dev authored and nyalldawson committed Nov 23, 2017
1 parent d36d4f8 commit 42bedc4f728985c12120b5989f8a3eafef7b79cd
@@ -815,8 +815,8 @@ void QgsRasterFileWriter::createVRT( int xSize, int ySize, const QgsCoordinateRe
if ( geoTransform )
{
QDomElement geoTransformElem = mVRTDocument.createElement( QStringLiteral( "GeoTransform" ) );
QString geoTransformString = QString::number( geoTransform[0] ) + ", " + QString::number( geoTransform[1] ) + ", " + QString::number( geoTransform[2] ) +
", " + QString::number( geoTransform[3] ) + ", " + QString::number( geoTransform[4] ) + ", " + QString::number( geoTransform[5] );
QString geoTransformString = QString::number( geoTransform[0], 'f', 6 ) + ", " + QString::number( geoTransform[1] ) + ", " + QString::number( geoTransform[2] ) +
", " + QString::number( geoTransform[3], 'f', 6 ) + ", " + QString::number( geoTransform[4] ) + ", " + QString::number( geoTransform[5] );
QDomText geoTransformText = mVRTDocument.createTextNode( geoTransformString );
geoTransformElem.appendChild( geoTransformText );
VRTDatasetElem.appendChild( geoTransformElem );
@@ -50,6 +50,7 @@ class TestQgsRasterFileWriter: public QObject
void writeTest();
void testCreateOneBandRaster();
void testCreateMultiBandRaster();
void testVrtCreation();
private:
bool writeTest( const QString &rasterName );
void log( const QString &msg );
@@ -271,6 +272,52 @@ void TestQgsRasterFileWriter::testCreateMultiBandRaster()
delete rlayer;
}

void TestQgsRasterFileWriter::testVrtCreation()
{
//create a raster layer that will be used in all tests...
QString srcFileName = mTestDataDir + QStringLiteral( "ALLINGES_RGF93_CC46_1_1.tif" );
QFileInfo rasterFileInfo( srcFileName );
std::unique_ptr< QgsRasterLayer > srcRasterLayer = qgis::make_unique< QgsRasterLayer >( rasterFileInfo.absoluteFilePath(), rasterFileInfo.completeBaseName() );

QTemporaryDir dir;
std::unique_ptr< QgsRasterFileWriter > rasterFileWriter = qgis::make_unique< QgsRasterFileWriter >( dir.path() + '/' + rasterFileInfo.completeBaseName() );

//2. Definition of the pyramid levels
QList<int> levelList;
levelList << 2 << 4 << 8 << 16 << 32 << 64 << 128;
rasterFileWriter->setPyramidsList( levelList );
//3. Pyramid format
rasterFileWriter->setPyramidsFormat( QgsRaster::PyramidsGTiff );
//4. Resampling method
rasterFileWriter->setPyramidsResampling( QStringLiteral( "NEAREST" ) );
//5. Tiled mode => true for vrt creation
rasterFileWriter->setTiledMode( true );
//6. Tile size
rasterFileWriter->setMaxTileWidth( 500 );
rasterFileWriter->setMaxTileHeight( 500 );
//7. Coordinate Reference System
QgsCoordinateReferenceSystem crs;
crs.createFromString( "EPSG:3946" );
//8. Prepare raster pipe
QgsRasterPipe pipe;
pipe.set( srcRasterLayer->dataProvider()->clone() );
// Let's do it !
QgsRasterFileWriter::WriterError res = rasterFileWriter->writeRaster( &pipe, srcRasterLayer->width(), srcRasterLayer->height(), srcRasterLayer->extent(), crs );
QCOMPARE( res, QgsRasterFileWriter::NoError );

// Now let's compare the georef of the original raster with the georef of the generated vrt file
std::unique_ptr< QgsRasterLayer > vrtRasterLayer = qgis::make_unique< QgsRasterLayer >( dir.path() + '/' + rasterFileInfo.completeBaseName() + '/' + rasterFileInfo.completeBaseName() + QStringLiteral( ".vrt" ), rasterFileInfo.completeBaseName() );

double xminVrt = vrtRasterLayer->extent().xMinimum();
double yminVrt = vrtRasterLayer->extent().yMaximum();
double xminOriginal = srcRasterLayer->extent().xMinimum();
double yminOriginal = srcRasterLayer->extent().yMaximum();

// Let's check if the georef of the original raster with the georef of the generated vrt file
QGSCOMPARENEAR( xminVrt, xminOriginal, srcRasterLayer->rasterUnitsPerPixelX() / 4 );
QGSCOMPARENEAR( yminVrt, yminOriginal, srcRasterLayer->rasterUnitsPerPixelY() / 4 );
}

void TestQgsRasterFileWriter::log( const QString &msg )
{
mReport += msg + "<br>";
@@ -0,0 +1,6 @@
0.020000000000
0
0
-0.020000000000
1964499.715750000207
5242972.348840000108
Binary file not shown.

0 comments on commit 42bedc4

Please sign in to comment.
You can’t perform that action at this time.