QGIS API Documentation 3.37.0-Master (fdefdf9c27f)
qgsalignraster.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalignraster.cpp
3 --------------------------------------
4 Date : June 2015
5 Copyright : (C) 2015 by Martin Dobias
6 Email : wonder dot sk at gmail dot com
7 ***************************************************************************
8 * *
9 * This program is free software; you can redistribute it and/or modify *
10 * it under the terms of the GNU General Public License as published by *
11 * the Free Software Foundation; either version 2 of the License, or *
12 * (at your option) any later version. *
13 * *
14 ***************************************************************************/
15
16#include "qgsalignraster.h"
17
18#include <gdalwarper.h>
19#include <ogr_srs_api.h>
20#include <cpl_conv.h>
21#include <limits>
22
23#include <QPair>
24#include <QString>
25
27#include "qgsrectangle.h"
28
29
30static double ceil_with_tolerance( double value )
31{
32 if ( std::fabs( value - std::round( value ) ) < 1e-6 )
33 return std::round( value );
34 else
35 return std::ceil( value );
36}
37
38static double floor_with_tolerance( double value )
39{
40 if ( std::fabs( value - std::round( value ) ) < 1e-6 )
41 return std::round( value );
42 else
43 return std::floor( value );
44}
45
46static double fmod_with_tolerance( double num, double denom )
47{
48 return num - floor_with_tolerance( num / denom ) * denom;
49}
50
51static QgsRectangle transform_to_extent( const double *geotransform, double xSize, double ySize )
52{
53 QgsRectangle r( geotransform[0],
54 geotransform[3],
55 geotransform[0] + geotransform[1] * xSize,
56 geotransform[3] + geotransform[5] * ySize );
57 r.normalize();
58 return r;
59}
60
61static int CPL_STDCALL _progress( double dfComplete, const char *pszMessage, void *pProgressArg )
62{
63 Q_UNUSED( pszMessage )
64
65 QgsAlignRaster::ProgressHandler *handler = ( ( QgsAlignRaster * ) pProgressArg )->progressHandler();
66 if ( handler )
67 return handler->progress( dfComplete );
68 else
69 return true;
70}
71
72static CPLErr rescalePreWarpChunkProcessor( void *pKern, void *pArg )
73{
74 GDALWarpKernel *kern = ( GDALWarpKernel * ) pKern;
75 const double cellsize = ( ( double * )pArg )[0];
76
77 for ( int nBand = 0; nBand < kern->nBands; ++nBand )
78 {
79 float *bandData = ( float * ) kern->papabySrcImage[nBand];
80 for ( int nLine = 0; nLine < kern->nSrcYSize; ++nLine )
81 {
82 for ( int nPixel = 0; nPixel < kern->nSrcXSize; ++nPixel )
83 {
84 bandData[nPixel] /= cellsize;
85 }
86 bandData += kern->nSrcXSize;
87 }
88 }
89 return CE_None;
90}
91
92static CPLErr rescalePostWarpChunkProcessor( void *pKern, void *pArg )
93{
94 GDALWarpKernel *kern = ( GDALWarpKernel * ) pKern;
95 const double cellsize = ( ( double * )pArg )[1];
96
97 for ( int nBand = 0; nBand < kern->nBands; ++nBand )
98 {
99 float *bandData = ( float * ) kern->papabyDstImage[nBand];
100 for ( int nLine = 0; nLine < kern->nDstYSize; ++nLine )
101 {
102 for ( int nPixel = 0; nPixel < kern->nDstXSize; ++nPixel )
103 {
104 bandData[nPixel] *= cellsize;
105 }
106 bandData += kern->nDstXSize;
107 }
108 }
109 return CE_None;
110}
111
112
113
115{
116 // parameters
119 mClipExtent[0] = mClipExtent[1] = mClipExtent[2] = mClipExtent[3] = 0;
120
121 // derived variables
122 mXSize = mYSize = 0;
123 for ( int i = 0; i < 6; ++i )
124 mGeoTransform[i] = 0;
125}
126
127void QgsAlignRaster::setClipExtent( double xmin, double ymin, double xmax, double ymax )
128{
129 mClipExtent[0] = xmin;
130 mClipExtent[1] = ymin;
131 mClipExtent[2] = xmax;
132 mClipExtent[3] = ymax;
133}
134
136{
137 setClipExtent( extent.xMinimum(), extent.yMinimum(),
138 extent.xMaximum(), extent.yMaximum() );
139}
140
142{
143 return QgsRectangle( mClipExtent[0], mClipExtent[1],
144 mClipExtent[2], mClipExtent[3] );
145}
146
147
148bool QgsAlignRaster::setParametersFromRaster( const QString &filename, const QString &destWkt, QSizeF customCellSize, QPointF customGridOffset )
149{
150 return setParametersFromRaster( RasterInfo( filename ), destWkt, customCellSize, customGridOffset );
151}
152
153bool QgsAlignRaster::setParametersFromRaster( const RasterInfo &rasterInfo, const QString &customCRSWkt, QSizeF customCellSize, QPointF customGridOffset )
154{
155 if ( customCRSWkt.isEmpty() || customCRSWkt == rasterInfo.crs() )
156 {
157 // use ref. layer to init input
158 mCrsWkt = rasterInfo.crs();
159
160 if ( !customCellSize.isValid() )
161 {
162 mCellSizeX = rasterInfo.cellSize().width();
163 mCellSizeY = rasterInfo.cellSize().height();
164 }
165 else
166 {
167 mCellSizeX = customCellSize.width();
168 mCellSizeY = customCellSize.height();
169 }
170
171 if ( customGridOffset.x() < 0 || customGridOffset.y() < 0 )
172 {
173 if ( !customCellSize.isValid() )
174 {
175 // using original raster's grid offset to be aligned with origin
176 mGridOffsetX = rasterInfo.gridOffset().x();
177 mGridOffsetY = rasterInfo.gridOffset().y();
178 }
179 else
180 {
181 // if using custom cell size: offset so that we are aligned
182 // with the original raster's origin point
183 mGridOffsetX = fmod_with_tolerance( rasterInfo.origin().x(), customCellSize.width() );
184 mGridOffsetY = fmod_with_tolerance( rasterInfo.origin().y(), customCellSize.height() );
185 }
186 }
187 else
188 {
189 mGridOffsetX = customGridOffset.x();
190 mGridOffsetY = customGridOffset.y();
191 }
192 }
193 else
194 {
195 QSizeF cs;
196 QPointF go;
197 if ( !suggestedWarpOutput( rasterInfo, customCRSWkt, &cs, &go ) )
198 {
199 mCrsWkt = QStringLiteral( "_error_" );
202 return false;
203 }
204
205 mCrsWkt = customCRSWkt;
206
207 if ( !customCellSize.isValid() )
208 {
209 mCellSizeX = cs.width();
210 mCellSizeY = cs.height();
211 }
212 else
213 {
214 mCellSizeX = customCellSize.width();
215 mCellSizeY = customCellSize.height();
216 }
217
218 if ( customGridOffset.x() < 0 || customGridOffset.y() < 0 )
219 {
220 mGridOffsetX = go.x();
221 mGridOffsetY = go.y();
222 }
223 else
224 {
225 mGridOffsetX = customGridOffset.x();
226 mGridOffsetY = customGridOffset.x();
227 }
228 }
229 return true;
230}
231
232
234{
235 mErrorMessage.clear();
236
237 if ( mCrsWkt == QLatin1String( "_error_" ) )
238 {
239 mErrorMessage = QObject::tr( "Unable to reproject." );
240 return false;
241 }
242
243 if ( mCellSizeX == 0 || mCellSizeY == 0 )
244 {
245 mErrorMessage = QObject::tr( "Cell size must not be zero." );
246 return false;
247 }
248
249 mXSize = mYSize = 0;
250 std::fill_n( mGeoTransform, 6, 0 );
251
252 double finalExtent[4] = { 0, 0, 0, 0 };
253
254 // for each raster: determine their extent in projected cfg
255 for ( int i = 0; i < mRasters.count(); ++i )
256 {
257 Item &r = mRasters[i];
258
259 const RasterInfo info( r.inputFilename );
260
261 QSizeF cs;
262 QgsRectangle extent;
263 if ( !suggestedWarpOutput( info, mCrsWkt, &cs, nullptr, &extent ) )
264 {
265 mErrorMessage = QString( "Failed to get suggested warp output.\n\n"
266 "File:\n%1\n\n"
267 "Source WKT:\n%2\n\nDestination WKT:\n%3" )
268 .arg( r.inputFilename,
269 info.mCrsWkt,
270 mCrsWkt );
271 return false;
272 }
273
274 r.srcCellSizeInDestCRS = cs.width() * cs.height();
275
276 if ( finalExtent[0] == 0 && finalExtent[1] == 0 && finalExtent[2] == 0 && finalExtent[3] == 0 )
277 {
278 // initialize with the first layer
279 finalExtent[0] = extent.xMinimum();
280 finalExtent[1] = extent.yMinimum();
281 finalExtent[2] = extent.xMaximum();
282 finalExtent[3] = extent.yMaximum();
283 }
284 else
285 {
286 // use intersection of rects
287 if ( extent.xMinimum() > finalExtent[0] ) finalExtent[0] = extent.xMinimum();
288 if ( extent.yMinimum() > finalExtent[1] ) finalExtent[1] = extent.yMinimum();
289 if ( extent.xMaximum() < finalExtent[2] ) finalExtent[2] = extent.xMaximum();
290 if ( extent.yMaximum() < finalExtent[3] ) finalExtent[3] = extent.yMaximum();
291 }
292 }
293
294 // count in extra clip extent (if present)
295 // 1. align requested rect to the grid - extend the rect if necessary
296 // 2. intersect with clip extent with final extent
297
298 if ( !( mClipExtent[0] == 0 && mClipExtent[1] == 0 && mClipExtent[2] == 0 && mClipExtent[3] == 0 ) )
299 {
300 // extend clip extent to grid
301 const double clipX0 = floor_with_tolerance( ( mClipExtent[0] - mGridOffsetX ) / mCellSizeX ) * mCellSizeX + mGridOffsetX;
302 const double clipY0 = floor_with_tolerance( ( mClipExtent[1] - mGridOffsetY ) / mCellSizeY ) * mCellSizeY + mGridOffsetY;
303 const double clipX1 = ceil_with_tolerance( ( mClipExtent[2] - clipX0 ) / mCellSizeX ) * mCellSizeX + clipX0;
304 const double clipY1 = ceil_with_tolerance( ( mClipExtent[3] - clipY0 ) / mCellSizeY ) * mCellSizeY + clipY0;
305 if ( clipX0 > finalExtent[0] ) finalExtent[0] = clipX0;
306 if ( clipY0 > finalExtent[1] ) finalExtent[1] = clipY0;
307 if ( clipX1 < finalExtent[2] ) finalExtent[2] = clipX1;
308 if ( clipY1 < finalExtent[3] ) finalExtent[3] = clipY1;
309 }
310
311 // align to grid - shrink the rect if necessary
312 // output raster grid configuration (with no rotation/shear)
313 // ... and raster width/height
314
315 const double originX = ceil_with_tolerance( ( finalExtent[0] - mGridOffsetX ) / mCellSizeX ) * mCellSizeX + mGridOffsetX;
316 const double originY = ceil_with_tolerance( ( finalExtent[1] - mGridOffsetY ) / mCellSizeY ) * mCellSizeY + mGridOffsetY;
317 const int xSize = floor_with_tolerance( ( finalExtent[2] - originX ) / mCellSizeX );
318 const int ySize = floor_with_tolerance( ( finalExtent[3] - originY ) / mCellSizeY );
319
320 if ( xSize <= 0 || ySize <= 0 )
321 {
322 mErrorMessage = QObject::tr( "No common intersecting area." );
323 return false;
324 }
325
326 mXSize = xSize;
327 mYSize = ySize;
328
329 // build final geotransform...
330 mGeoTransform[0] = originX;
332 mGeoTransform[2] = 0;
333 mGeoTransform[3] = originY + ( mCellSizeY * ySize );
334 mGeoTransform[4] = 0;
336
337 return true;
338}
339
340
342{
343 return QSize( mXSize, mYSize );
344}
345
347{
348 return transform_to_extent( mGeoTransform, mXSize, mYSize );
349}
350
351
353{
354 mErrorMessage.clear();
355
356 // consider extent of all layers and setup geotransform and output grid size
357 if ( !checkInputParameters() )
358 return false;
359
360 //dump();
361
362 const auto constMRasters = mRasters;
363 for ( const Item &r : constMRasters )
364 {
365 if ( !createAndWarp( r ) )
366 return false;
367 }
368 return true;
369}
370
371
373{
374 qDebug( "---ALIGN------------------" );
375 qDebug( "wkt %s", mCrsWkt.toLatin1().constData() );
376 qDebug( "w/h %d,%d", mXSize, mYSize );
377 qDebug( "transform" );
378 qDebug( "%6.2f %6.2f %6.2f", mGeoTransform[0], mGeoTransform[1], mGeoTransform[2] );
379 qDebug( "%6.2f %6.2f %6.2f", mGeoTransform[3], mGeoTransform[4], mGeoTransform[5] );
380
381 const QgsRectangle e = transform_to_extent( mGeoTransform, mXSize, mYSize );
382 qDebug( "extent %s", e.toString().toLatin1().constData() );
383}
384
386{
387 int bestIndex = -1;
388 double bestCellArea = INFINITY;
389 QSizeF cs;
390 int i = 0;
391
392 // using WGS84 as a destination CRS... but maybe some projected CRS
393 // would be a better a choice to more accurately compute areas?
394 // (Why earth is not flat???)
395 const QgsCoordinateReferenceSystem destCRS( QStringLiteral( "EPSG:4326" ) );
396 const QString destWkt = destCRS.toWkt( Qgis::CrsWktVariant::PreferredGdal );
397
398 const auto constMRasters = mRasters;
399 for ( const Item &raster : constMRasters )
400 {
401 if ( !suggestedWarpOutput( RasterInfo( raster.inputFilename ), destWkt, &cs ) )
402 return false;
403
404 const double cellArea = cs.width() * cs.height();
405 if ( cellArea < bestCellArea )
406 {
407 bestCellArea = cellArea;
408 bestIndex = i;
409 }
410 ++i;
411 }
412
413 return bestIndex;
414}
415
416
418{
419 GDALDriverH hDriver = GDALGetDriverByName( "GTiff" );
420 if ( !hDriver )
421 {
422 mErrorMessage = QStringLiteral( "GDALGetDriverByName(GTiff) failed." );
423 return false;
424 }
425
426 // Open the source file.
427 const gdal::dataset_unique_ptr hSrcDS( GDALOpen( raster.inputFilename.toUtf8().constData(), GA_ReadOnly ) );
428 if ( !hSrcDS )
429 {
430 mErrorMessage = QObject::tr( "Unable to open input file: %1" ).arg( raster.inputFilename );
431 return false;
432 }
433
434 // Create output with same datatype as first input band.
435
436 const int bandCount = GDALGetRasterCount( hSrcDS.get() );
437 const GDALDataType eDT = GDALGetRasterDataType( GDALGetRasterBand( hSrcDS.get(), 1 ) );
438
439 // Create the output file.
440 const gdal::dataset_unique_ptr hDstDS( GDALCreate( hDriver, raster.outputFilename.toUtf8().constData(), mXSize, mYSize,
441 bandCount, eDT, nullptr ) );
442 if ( !hDstDS )
443 {
444 mErrorMessage = QObject::tr( "Unable to create output file: %1" ).arg( raster.outputFilename );
445 return false;
446 }
447
448 // Write out the projection definition.
449 GDALSetProjection( hDstDS.get(), mCrsWkt.toLatin1().constData() );
450 GDALSetGeoTransform( hDstDS.get(), mGeoTransform );
451
452 // Copy the color table, if required.
453 GDALColorTableH hCT = GDALGetRasterColorTable( GDALGetRasterBand( hSrcDS.get(), 1 ) );
454 if ( hCT )
455 GDALSetRasterColorTable( GDALGetRasterBand( hDstDS.get(), 1 ), hCT );
456
457 // -----------------------------------------------------------------------
458
459 // Setup warp options.
460 gdal::warp_options_unique_ptr psWarpOptions( GDALCreateWarpOptions() );
461 psWarpOptions->hSrcDS = hSrcDS.get();
462 psWarpOptions->hDstDS = hDstDS.get();
463
464 psWarpOptions->nBandCount = GDALGetRasterCount( hSrcDS.get() );
465 psWarpOptions->panSrcBands = ( int * ) CPLMalloc( sizeof( int ) * psWarpOptions->nBandCount );
466 psWarpOptions->panDstBands = ( int * ) CPLMalloc( sizeof( int ) * psWarpOptions->nBandCount );
467 for ( int i = 0; i < psWarpOptions->nBandCount; ++i )
468 {
469 psWarpOptions->panSrcBands[i] = i + 1;
470 psWarpOptions->panDstBands[i] = i + 1;
471 }
472
473 psWarpOptions->eResampleAlg = static_cast< GDALResampleAlg >( raster.resampleMethod );
474
475 // our progress function
476 psWarpOptions->pfnProgress = _progress;
477 psWarpOptions->pProgressArg = this;
478
479 // Establish reprojection transformer.
480 psWarpOptions->pTransformerArg =
481 GDALCreateGenImgProjTransformer( hSrcDS.get(), GDALGetProjectionRef( hSrcDS.get() ),
482 hDstDS.get(), GDALGetProjectionRef( hDstDS.get() ),
483 FALSE, 0.0, 1 );
484 psWarpOptions->pfnTransformer = GDALGenImgProjTransform;
485
486 double rescaleArg[2];
487 if ( raster.rescaleValues )
488 {
489 rescaleArg[0] = raster.srcCellSizeInDestCRS; // source cell size
490 rescaleArg[1] = mCellSizeX * mCellSizeY; // destination cell size
491 psWarpOptions->pfnPreWarpChunkProcessor = rescalePreWarpChunkProcessor;
492 psWarpOptions->pfnPostWarpChunkProcessor = rescalePostWarpChunkProcessor;
493 psWarpOptions->pPreWarpProcessorArg = rescaleArg;
494 psWarpOptions->pPostWarpProcessorArg = rescaleArg;
495 // force use of float32 data type as that is what our pre/post-processor uses
496 psWarpOptions->eWorkingDataType = GDT_Float32;
497 }
498
499 // Initialize and execute the warp operation.
500 GDALWarpOperation oOperation;
501 oOperation.Initialize( psWarpOptions.get() );
502 oOperation.ChunkAndWarpImage( 0, 0, mXSize, mYSize );
503
504 GDALDestroyGenImgProjTransformer( psWarpOptions->pTransformerArg );
505 return true;
506}
507
508bool QgsAlignRaster::suggestedWarpOutput( const QgsAlignRaster::RasterInfo &info, const QString &destWkt, QSizeF *cellSize, QPointF *gridOffset, QgsRectangle *rect )
509{
510 // Create a transformer that maps from source pixel/line coordinates
511 // to destination georeferenced coordinates (not destination
512 // pixel line). We do that by omitting the destination dataset
513 // handle (setting it to nullptr).
514 void *hTransformArg = GDALCreateGenImgProjTransformer( info.mDataset.get(), info.mCrsWkt.toLatin1().constData(), nullptr, destWkt.toLatin1().constData(), FALSE, 0, 1 );
515 if ( !hTransformArg )
516 return false;
517
518 // Get approximate output georeferenced bounds and resolution for file.
519 double adfDstGeoTransform[6];
520 double extents[4];
521 int nPixels = 0, nLines = 0;
522 CPLErr eErr;
523 eErr = GDALSuggestedWarpOutput2( info.mDataset.get(),
524 GDALGenImgProjTransform, hTransformArg,
525 adfDstGeoTransform, &nPixels, &nLines, extents, 0 );
526 GDALDestroyGenImgProjTransformer( hTransformArg );
527
528 if ( eErr != CE_None )
529 return false;
530
531 const QSizeF cs( std::fabs( adfDstGeoTransform[1] ), std::fabs( adfDstGeoTransform[5] ) );
532
533 if ( rect )
534 *rect = QgsRectangle( extents[0], extents[1], extents[2], extents[3] );
535 if ( cellSize )
536 *cellSize = cs;
537 if ( gridOffset )
538 *gridOffset = QPointF( fmod_with_tolerance( adfDstGeoTransform[0], cs.width() ),
539 fmod_with_tolerance( adfDstGeoTransform[3], cs.height() ) );
540 return true;
541}
542
543
544//----------
545
546
547QgsAlignRaster::RasterInfo::RasterInfo( const QString &layerpath )
548{
549 mDataset.reset( GDALOpen( layerpath.toUtf8().constData(), GA_ReadOnly ) );
550 if ( !mDataset )
551 return;
552
553 mXSize = GDALGetRasterXSize( mDataset.get() );
554 mYSize = GDALGetRasterYSize( mDataset.get() );
555
556 ( void ) GDALGetGeoTransform( mDataset.get(), mGeoTransform );
557
558 // TODO: may be null or empty string
559 mCrsWkt = QString::fromLatin1( GDALGetProjectionRef( mDataset.get() ) );
560
561 mBandCnt = GDALGetBandNumber( mDataset.get() );
562}
563
565{
566 return QSizeF( std::fabs( mGeoTransform[1] ), std::fabs( mGeoTransform[5] ) );
567}
568
570{
571 return QPointF( fmod_with_tolerance( mGeoTransform[0], cellSize().width() ),
572 fmod_with_tolerance( mGeoTransform[3], cellSize().height() ) );
573}
574
576{
577 return transform_to_extent( mGeoTransform, mXSize, mYSize );
578}
579
581{
582 return QPointF( mGeoTransform[0], mGeoTransform[3] );
583}
584
586{
587 qDebug( "---RASTER INFO------------------" );
588 qDebug( "wkt %s", mCrsWkt.toLatin1().constData() );
589 qDebug( "w/h %d,%d", mXSize, mYSize );
590 qDebug( "cell x/y %f,%f", cellSize().width(), cellSize().width() );
591
592 const QgsRectangle r = extent();
593 qDebug( "extent %s", r.toString().toLatin1().constData() );
594
595 qDebug( "transform" );
596 qDebug( "%6.2f %6.2f %6.2f", mGeoTransform[0], mGeoTransform[1], mGeoTransform[2] );
597 qDebug( "%6.2f %6.2f %6.2f", mGeoTransform[3], mGeoTransform[4], mGeoTransform[5] );
598}
599
600double QgsAlignRaster::RasterInfo::identify( double mx, double my )
601{
602 GDALRasterBandH hBand = GDALGetRasterBand( mDataset.get(), 1 );
603
604 // must not be rotated in order for this to work
605 const int px = int( ( mx - mGeoTransform[0] ) / mGeoTransform[1] );
606 const int py = int( ( my - mGeoTransform[3] ) / mGeoTransform[5] );
607
608 float *pafScanline = ( float * ) CPLMalloc( sizeof( float ) );
609 const CPLErr err = GDALRasterIO( hBand, GF_Read, px, py, 1, 1,
610 pafScanline, 1, 1, GDT_Float32, 0, 0 );
611 const double value = err == CE_None ? pafScanline[0] : std::numeric_limits<double>::quiet_NaN();
612 CPLFree( pafScanline );
613
614 return value;
615}
616
@ PreferredGdal
Preferred format for conversion of CRS to WKT for use with the GDAL library.
QgsAlignRaster takes one or more raster layers and warps (resamples) them so they have the same:
void dump() const
write contents of the object to standard error stream - for debugging
int mYSize
Computed raster grid height.
bool setParametersFromRaster(const RasterInfo &rasterInfo, const QString &customCRSWkt=QString(), QSizeF customCellSize=QSizeF(), QPointF customGridOffset=QPointF(-1, -1))
Set destination CRS, cell size and grid offset from a raster file.
QSizeF cellSize() const
Gets output cell size.
QgsRectangle clipExtent() const
Gets clipping extent (region of interest).
QPointF gridOffset() const
double mGridOffsetX
Destination grid offset - expected to be in interval <0,cellsize)
QgsRectangle alignedRasterExtent() const
Returns the expected extent of the resulting aligned raster.
bool run()
Run the alignment process.
double mClipExtent[4]
Optional clip extent: sets "requested area" which be extended to fit the raster grid.
List mRasters
List of rasters to be aligned (with their output files and other options)
bool createAndWarp(const Item &raster)
Internal function for processing of one raster (1. create output, 2. do the alignment)
QString mCrsWkt
Destination CRS - stored in well-known text (WKT) format.
void setClipExtent(double xmin, double ymin, double xmax, double ymax)
Configure clipping extent (region of interest).
int mXSize
Computed raster grid width.
int suggestedReferenceLayer() const
Returns the index of the layer which has smallest cell size (returns -1 on error)
bool checkInputParameters()
Determine destination extent from the input rasters and calculate derived values.
double mCellSizeX
Destination cell size.
double mGeoTransform[6]
Computed geo-transform.
static bool suggestedWarpOutput(const RasterInfo &info, const QString &destWkt, QSizeF *cellSize=nullptr, QPointF *gridOffset=nullptr, QgsRectangle *rect=nullptr)
Determine suggested output of raster warp to a different CRS. Returns true on success.
QSize alignedRasterSize() const
Returns the expected size of the resulting aligned raster.
QString mErrorMessage
Last error message from run()
This class represents a coordinate reference system (CRS).
QString toWkt(Qgis::CrsWktVariant variant=Qgis::CrsWktVariant::Wkt1Gdal, bool multiline=false, int indentationWidth=4) const
Returns a WKT representation of this CRS.
A rectangle specified with double values.
Definition: qgsrectangle.h:42
QString toString(int precision=16) const
Returns a string representation of form xmin,ymin : xmax,ymax Coordinates will be truncated to the sp...
double xMinimum() const
Returns the x minimum value (left side of rectangle).
Definition: qgsrectangle.h:201
double yMinimum() const
Returns the y minimum value (bottom side of rectangle).
Definition: qgsrectangle.h:211
double xMaximum() const
Returns the x maximum value (right side of rectangle).
Definition: qgsrectangle.h:196
double yMaximum() const
Returns the y maximum value (top side of rectangle).
Definition: qgsrectangle.h:206
std::unique_ptr< std::remove_pointer< GDALDatasetH >::type, GDALDatasetCloser > dataset_unique_ptr
Scoped GDAL dataset.
Definition: qgsogrutils.h:157
std::unique_ptr< GDALWarpOptions, GDALWarpOptionsDeleter > warp_options_unique_ptr
Scoped GDAL warp options.
Definition: qgsogrutils.h:172
Definition of one raster layer for alignment.
Qgis::GdalResampleAlgorithm resampleMethod
resampling method to be used
bool rescaleValues
rescaling of values according to the change of pixel size
double srcCellSizeInDestCRS
used for rescaling of values (if necessary)
QString inputFilename
filename of the source raster
QString outputFilename
filename of the newly created aligned raster (will be overwritten if exists already)
Helper struct to be sub-classed for progress reporting.
virtual bool progress(double complete)=0
Method to be overridden for progress reporting.
Utility class for gathering information about rasters.
QString crs() const
Returns the CRS in WKT format.
int mYSize
raster grid size Y
gdal::dataset_unique_ptr mDataset
handle to open GDAL dataset
void dump() const
Write contents of the object to standard error stream - for debugging.
RasterInfo(const QString &layerpath)
Construct raster info with a path to a raster file.
double mGeoTransform[6]
geotransform coefficients
QgsRectangle extent() const
Returns the extent of the raster.
QPointF gridOffset() const
Returns the grid offset.
QSizeF cellSize() const
Returns the cell size in map units.
double identify(double mx, double my)
Gets raster value at the given coordinates (from the first band)
int mXSize
raster grid size X
QPointF origin() const
Returns the origin of the raster.
int mBandCnt
number of raster's bands
QString mCrsWkt
CRS stored in WKT format.