qgis / QGIS

QGIS is a free, open source, cross platform (lin/win/mac) geographical information system (GIS)
https://qgis.org
GNU General Public License v2.0
10.11k stars 2.93k forks source link

GDAL Raster calculator giving incorrect results when confronted with certain Raster layers #42871

Closed iiepdev closed 3 years ago

iiepdev commented 3 years ago

Describe the bug

I am using GDAL Raster calculator to automatically reclassify Raster layers so that regardless of their initial values, they are transformed to be between 0 and 1 (or 0 and 100, as needed).

Since this needs to happen within a Graphical model, the algorithm needs to be able to automatically detect the minimum and maximum values within each Raster, and then make the following transformation:

NewRasterValue = (RasterValue - min(RasterValue))/(max(RasterValue) - min(RasterValue))

This will automatically shift the distribution to be bounded by 0 and 1 (and can be then multiplied by 100 to be scaled to be 0 to 100).

To achieve this I'm using the GDAL Raster calculator, with the following parameters:

Input layer A: Raster file Number of raster band for A: Band 1 Calculation in gdalnumeric syntax using +-/ or any numpy array functions: (A - amin(A))/(amax(A) - amin(A))100 Output raster type: Float32

When I do this for some raster files, I get the expected results. However, with some other layers the GDAL Raster calculator outputs a raster where every pixel is equal to 1 (or to a 100, in the case of the equation above).

How to Reproduce

  1. Raster files to use: This one doesn't result in the error: https://box.iiep.unesco.org/s/MgHkfS3bYyp8QWr This one results in the error: https://box.iiep.unesco.org/s/NCRmKj4oPSc7DBo
  2. Select GDAL's Raster calculator on the Processing Toolbox
  3. Apply the parameters presented in the previous section
  4. See error --> With one of the raster files, the results are as expected. With the other, the value of all pixels becomes 100.

QGIS and OS versions

QGIS version 3.16.0-Hannover QGIS code revision 43b64b13f3
Compiled against Qt 5.11.2 Running against Qt 5.11.2
Compiled against GDAL/OGR 3.1.4 Running against GDAL/OGR 3.1.4
Compiled against GEOS 3.8.1-CAPI-1.13.3 Running against GEOS 3.8.1-CAPI-1.13.3
Compiled against SQLite 3.29.0 Running against SQLite 3.29.0
PostgreSQL Client Version 11.5 SpatiaLite Version 4.3.0
QWT Version 6.1.3 QScintilla2 Version 2.10.8
Compiled against PROJ 6.3.2 Running against PROJ Rel. 6.3.2, May 1st, 2020
OS Version Windows 10 (10.0)
Active python plugins BivariateLegend; DataPlotly; geetimeseriesexplorer; GlobeBuilder; minimum_spanning_tree; mmqgis; ORStools; pluginbuilder3; plugin_reloader; processing_r; qgis2web; QNEAT3; QuickOSM; sprague_multipliers; db_manager; processing

Additional context

gioman commented 3 years ago

GDAL Raster calculator

@iiepdev you must try the same operation/command from the command line: do you get the same results? be sure you are using the same GDAL version/installation being used by QGIS.

iiepdev commented 3 years ago

GDAL Raster calculator

@iiepdev you must try the same operation/command from the command line: do you get the same results? be sure you are using the same GDAL version/installation being used by QGIS.

Thank you @gioman

I have never used GDAL from the command line. I imagine you mean to run this command from the OSGeo4W Shell or from the Command prompt, using the GDAL/OGR console call command (which for this case is: gdal_calc --calc "(A - amin(A))/(amax(A) - amin(A))*100" --format GTiff --type Float32 -A "C:/Users/g.vargas/BOX/IIEP_MyProjects/MP_01000298_RND_SDA/WorkFiles_Experts/298-Issue-Papers/298-Issue-Paper-UNOSAT/Geospatial Data - Part A/Part A/Raster Data/Hazard_Volcano_Aceh.tif" --A_band 1 --outfile C:/Users/g.vargas/AppData/Local/Temp/processing_MqvAmr/0713417ed361418f9af6f1e4d3f3ef6c/OUTPUT.tif).

When I try to run either option it doesn't recognize the command, which means I might need to install the python bindings for GDAL? (https://stackoverflow.com/questions/61901761/running-gdal-from-the-command-line-on-windows-10-x64)

Am I missing something?

gioman commented 3 years ago

@iiepdev yes the command is the one you refer, you must change it a little to be able to run it from the CLI manually, example:

gdal_calc --calc "(A - amin(A))/(amax(A) - amin(A))*100" --format GTiff --type Float32 -A "C:/Users/g.vargas/BOX/IIEP_MyProjects/MP_01000298_RND_SDA/WorkFiles_Experts/298-Issue-Papers/298-Issue-Paper-UNOSAT/Geospatial Data - Part A/Part A/Raster Data/Hazard_Volcano_Aceh.tif" --A_band 1 --outfile C:/Users/g.vargas/BOX/IIEP_MyProjects/MP_01000298_RND_SDA/WorkFiles_Experts/298-Issue-Papers/298-Issue-Paper-UNOSAT/Geospatial Data - Part A/Part A/OUTPUT.tif

When I try to run either option it doesn't recognize the command, which means I might need to install the python bindings for GDAL? (https://stackoverflow.com/questions/61901761/running-gdal-from-the-command-line-on-windows-10-x64)

no, just open the "ogseo4w" shell and run the command from it. That shell shortcut is inside the "QGIS" group that the QGIS installer created in the start menu of your Windows.

github-actions[bot] commented 3 years ago

The QGIS project highly values your report and would love to see it addressed. However, this issue has been left in feedback mode for the last 14 days and is being automatically marked as "stale". If you would like to continue with this issue, please provide any missing information or answer any open questions. If you could resolve the issue yourself meanwhile, please leave a note for future readers with the same problem and close the issue. In case you should have any uncertainty, please leave a comment and we will be happy to help you proceed with this issue. If there is no further activity on this issue, it will be closed in a week.

github-actions[bot] commented 3 years ago

While we hate to see this happen, this issue has been automatically closed because it has not had any activity in the last 42 days despite being marked as feedback. If this issue should be reconsidered, please follow the guidelines in the previous comment and reopen this issue. Or, if you have any further questions, there are also further support channels that can help you.