zoom to native resolution is not working correctly

[Andreas Janz]

Start the EnMAP-Box like so:

‌

from qgis._core import QgsRasterLayer

from enmapbox import EnMAPBox, initAll
from enmapbox.exampledata import enmap
from enmapbox.testing import start_app

qgsApp = start_app()
initAll()
enmapBox = EnMAPBox(None)

layer = QgsRasterLayer(enmap)
mapDock = enmapBox.onDataDropped([layer])
mapCanvas = mapDock.mapCanvas()
def debug():
    print(mapCanvas.mapUnitsPerPixel())
mapCanvas.extentsChanged.connect(debug)

qgsApp.exec_()

In the GUI, zoom to native resolution and have a look at the debug printouts. Instead of the expected 30 m, I get

30.303030303030305

This is almost correct, but when zooming to native resolution, I really want to see each image pixel on exactly one screen pixel. A difference of 0.3 m will introduce resampled pixels.

For all the old ENVI Classic folks (which includes me), this is critical.

[Andreas Janz]

resolves #853

→ <<cset 02db23914397>>

[Andreas Janz]

resolves #853

→ <<cset f8ce540458d8>>

[Andreas Janz]

resolves #853

→ <<cset fa5f1bd253c5>>

[Andreas Janz]

resolves #853

→ <<cset b1a8c201c652>>

[Benjamin Jakimow]

How should it work in case the canvas CRS and that of the raster layer are different?

Let the layer CRS be EPSG 32632 (UTM North, zone 32)
case1: canvas CRS is EPSG 32633 --> layer pixel grid slightly turned against the canvas CRS grid
case2: canvas EPSG 4326 (lat/lon) --> same case, depending on the spatial extent?

[Andreas Janz]

Have you tried it? Does is show wrong results?

Please provide a testdataset to reproduce the problem and a screenshot showing the wrong results.

I’m not sure if there really is a problem with the current implementation.

[Benjamin Jakimow]

addresses #853 (QPS updates) - added test_PixelScaleExtentMapTool - re-implemented zoom-to-native resolution from QGIS\src\app\qgisapp.cpp in maptools.py

→ <<cset d397c8c3ac36>>

[Benjamin Jakimow]

added test case with different pixel resolutions resolves #853

→ <<cset b5626658f09d>>

[Benjamin Jakimow]

I haven’t and it might work indeed also in strange combinations of layer and canvas CRS. However, I was wondering about the constant scale values and checked the QGIS implementation. It does not use any of them but refers to the spatial extent of the central map layer pixel instead (so I re-implemented, see b562665).

In particular for geographic canvas CRS I’d expect that mapCanvas.mapUnitsPerPixel() increases with higher latitudes.

void QgisApp::legendLayerZoomNative()
{
  if ( !mLayerTreeView )
    return;

  //find current Layer
  QgsMapLayer *currentLayer = mLayerTreeView->currentLayer();
  if ( !currentLayer )
    return;

  if ( QgsRasterLayer *layer = qobject_cast<QgsRasterLayer *>( currentLayer ) )
  {
    QgsDebugMsgLevel( "Raster units per pixel  : " + QString::number( layer->rasterUnitsPerPixelX() ), 2 );
    QgsDebugMsgLevel( "MapUnitsPerPixel before : " + QString::number( mMapCanvas->mapUnitsPerPixel() ), 2 );

    QList< double >nativeResolutions;
    if ( layer->dataProvider() )
    {
      nativeResolutions = layer->dataProvider()->nativeResolutions();
    }

    // get length of central canvas pixel width in source raster crs
    QgsRectangle e = mMapCanvas->extent();
    QSize s = mMapCanvas->mapSettings().outputSize();
    QgsPointXY p1( e.center().x(), e.center().y() );
    QgsPointXY p2( e.center().x() + e.width() / s.width(), e.center().y() + e.height() / s.height() );
    QgsCoordinateTransform ct( mMapCanvas->mapSettings().destinationCrs(), layer->crs(), QgsProject::instance() );
    p1 = ct.transform( p1 );
    p2 = ct.transform( p2 );
    const double diagonalSize = std::sqrt( p1.sqrDist( p2 ) ); // width (actually the diagonal) of reprojected pixel
    if ( !nativeResolutions.empty() )
    {
      // find closest native resolution
      QList< double > diagonalNativeResolutions;
      diagonalNativeResolutions.reserve( nativeResolutions.size() );
      for ( double d : nativeResolutions )
        diagonalNativeResolutions << std::sqrt( 2 * d * d );

      int i;
      for ( i = 0; i < diagonalNativeResolutions.size() && diagonalNativeResolutions.at( i ) < diagonalSize; i++ )
      {
        QgsDebugMsgLevel( QStringLiteral( "test resolution %1: %2" ).arg( i ).arg( diagonalNativeResolutions.at( i ) ), 2 );
      }
      if ( i == nativeResolutions.size() ||
           ( i > 0 && ( ( diagonalNativeResolutions.at( i ) - diagonalSize ) > ( diagonalSize - diagonalNativeResolutions.at( i - 1 ) ) ) ) )
      {
        QgsDebugMsgLevel( QStringLiteral( "previous resolution" ), 2 );
        i--;
      }

      mMapCanvas->zoomByFactor( nativeResolutions.at( i ) / mMapCanvas->mapUnitsPerPixel() );
    }
    else
    {
      mMapCanvas->zoomByFactor( std::sqrt( layer->rasterUnitsPerPixelX() * layer->rasterUnitsPerPixelX() + layer->rasterUnitsPerPixelY() * layer->rasterUnitsPerPixelY() ) / diagonalSize );
    }

    mMapCanvas->refresh();
    QgsDebugMsgLevel( "MapUnitsPerPixel after  : " + QString::number( mMapCanvas->mapUnitsPerPixel() ), 2 );
  }
}

‌

[Andreas Janz]

In particular for geographic canvas CRS I’d expect that mapCanvas.mapUnitsPerPixel() increases with higher latitudes.

If I remember correctly, QGIS (and also other software like GoogleEarthEngine) uses the size at the equator.

[Andreas Janz]

However, I was wondering about the constant scale values and checked the QGIS implementation. It does not use any of them but refers to the spatial extent of the central map layer pixel instead (so I re-implemented, see b562665).

Ah nice, I wasn’t able to find the C++ code for that

[Benjamin Jakimow]

Squashed 'enmapbox/qgispluginsupport/' changes from 697d895ef..8604ddec4

8604ddec4 flaked files f04bccfaf spectralprocessingdialog.py writes & reads spectral setting to/from temporary files 94a92b761 spectralprocessingdialog.py write spectral setting to temporary files 95fbb2ea1 Merge branch 'write_spectral_setting' 0171ee8da python-tests.yml added releast-3_22 to test matrix cf607f050 added workaround for missing QgsFieldCalculator 9b2199131 added workaround QgsProcessingParametersGenerator.Flags() 2d9c72c6b hides QgsFieldCalculator if not available 8886405dd added repaint trigger added typechecking benchmark a60f4d101 added repaint trigger added typechecking benchmark 56afe2b3d enhances SpectralSettings added release-3_22 to test routines Signed-off-by: jakimowb [e-mail] dccdef480 Spectral Processing Dialog got static window title, as algo name is shown in separate widget fixed XML read error cb6eb54b5 added QgsFieldCalculator 3b0199fd4 refactoring spectral profile plots included spectral processing dialog 03ab35255 WIP refactoring profile plot model towards using QStandardItems 0f6b2f813 addresses #1 - added SelectMapLayerDialog b1761414c addresses #1 - added LayerRendererVisualization + unit test - flaked repo 59e6a783f addresses #b902032 addresses Issue #875 improved spectral profile plotting added workaround for wavelength definition at provider level 121dcec80 added test case with different pixel resolutions resolves #853 ddd03c777 addresses #853 (QPS updates) - added test_PixelScaleExtentMapTool - re-implemented zoom-to-native resolution from QGIS\src\app\qgisapp.cpp in maptools.py f74e1a3da qgisinfo.py works without start_app speclibrasterdataprovider initialized during init.py phase 4259415eb addresses #905 : temporary profile candidate gets color from profile source panel a236b916f qps pyqtgraph - try-except itemChange

git-subtree-dir: enmapbox/qgispluginsupport git-subtree-split: 8604ddec47b6d4add389005ec6bc3a09a8bde3fd

→ <<cset 03db28f8ba8e>>

[Benjamin Jakimow]

Squashed 'enmapbox/qgispluginsupport/' changes from 697d895ef..8604ddec4

8604ddec4 flaked files f04bccfaf spectralprocessingdialog.py writes & reads spectral setting to/from temporary files 94a92b761 spectralprocessingdialog.py write spectral setting to temporary files 95fbb2ea1 Merge branch 'write_spectral_setting' 0171ee8da python-tests.yml added releast-3_22 to test matrix cf607f050 added workaround for missing QgsFieldCalculator 9b2199131 added workaround QgsProcessingParametersGenerator.Flags() 2d9c72c6b hides QgsFieldCalculator if not available 8886405dd added repaint trigger added typechecking benchmark a60f4d101 added repaint trigger added typechecking benchmark 56afe2b3d enhances SpectralSettings added release-3_22 to test routines Signed-off-by: jakimowb [e-mail] dccdef480 Spectral Processing Dialog got static window title, as algo name is shown in separate widget fixed XML read error cb6eb54b5 added QgsFieldCalculator 3b0199fd4 refactoring spectral profile plots included spectral processing dialog 03ab35255 WIP refactoring profile plot model towards using QStandardItems 0f6b2f813 addresses #1 - added SelectMapLayerDialog b1761414c addresses #1 - added LayerRendererVisualization + unit test - flaked repo 59e6a783f addresses #b902032 addresses Issue #875 improved spectral profile plotting added workaround for wavelength definition at provider level 121dcec80 added test case with different pixel resolutions resolves #853 ddd03c777 addresses #853 (QPS updates) - added test_PixelScaleExtentMapTool - re-implemented zoom-to-native resolution from QGIS\src\app\qgisapp.cpp in maptools.py f74e1a3da qgisinfo.py works without start_app speclibrasterdataprovider initialized during init.py phase 4259415eb addresses #905 : temporary profile candidate gets color from profile source panel a236b916f qps pyqtgraph - try-except itemChange

git-subtree-dir: enmapbox/qgispluginsupport git-subtree-split: 8604ddec47b6d4add389005ec6bc3a09a8bde3fd

→ <<cset 03db28f8ba8e>>