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