search

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.56k stars 3k forks source link

QGIS treats AAIGrid file as Float 64 even if I explicitly tell it not to #46798

Open AlisterH opened 2 years ago

AlisterH commented 2 years ago

What is the bug or the crash?

When loading an AAIGrid raster QGIS treats it as Float 64 even if I explicitly tell it not to.

If I "layer>Save as" to Geotiff, I get a Float64 file.

To save to Float32, the workaround is simply to use the gdal_translate or gdalwarp processing algorithms, which produce Float32 output even if I don't specify the output data type. (As far as I know there is no way to specify the output data type when using the "Save as" dialog - unless there's a creation option that could be used which I missed?)

Steps to reproduce the issue

Example file.zip

Unzip and drag and drop into QGIS and check the Layer properties, and you can see QGIS treats this file as Float64.

Instead of drag and dropping, use the data source manager, select the file and under Options specify the DATATYPE as Int32. Check the file properties, and again QGIS is treating it as Float64. Why do we have an option to specify DATATYPE if it is just ignored?

Gdal calls the file Float32, and treats it as such:

C:\Processing>gdalinfo test.asc
Driver: AAIGrid/Arc/Info ASCII Grid
Files: test.asc
Size is 95, 56
Origin = (415630.000000000000000,759380.000000000000000)
Pixel Size = (10.000000000000000,-10.000000000000000)
Corner Coordinates:
Upper Left  (  415630.000,  759380.000)
Lower Left  (  415630.000,  758820.000)
Upper Right (  416580.000,  759380.000)
Lower Right (  416580.000,  758820.000)
Center      (  416105.000,  759100.000)
Band 1 Block=95x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
  Metadata:
    STATISTICS_MAXIMUM=78.616
    STATISTICS_MEAN=71.720051459052
    STATISTICS_MINIMUM=45.724
    STATISTICS_STDDEV=5.7599521748039
    STATISTICS_VALID_PERCENT=59.91

C:\Processing>gdal_translate test.asc test.tif
Input file size is 95, 56
0...10...20...30...40...50...60...70...80...90...100 - done.

C:\Processing>gdalinfo test.tif
Driver: GTiff/GeoTIFF
Files: test.tif
Size is 95, 56
Origin = (415630.000000000000000,759380.000000000000000)
Pixel Size = (10.000000000000000,-10.000000000000000)
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  415630.000,  759380.000)
Lower Left  (  415630.000,  758820.000)
Upper Right (  416580.000,  759380.000)
Lower Right (  416580.000,  758820.000)
Center      (  416105.000,  759100.000)
Band 1 Block=95x21 Type=Float32, ColorInterp=Gray
  Min=45.724 Max=78.616
  Minimum=45.724, Maximum=78.616, Mean=71.720, StdDev=5.760
  NoData Value=-9999
  Metadata:
    STATISTICS_MAXIMUM=78.616
    STATISTICS_MEAN=71.720051459052
    STATISTICS_MINIMUM=45.724
    STATISTICS_STDDEV=5.7599521748039
    STATISTICS_VALID_PERCENT=59.91

Obviously it isn't really either Float32 or Float64, but if we want to do something like save as to TIFF, Float32 is more than accurate enough, and Float64 just bloats the file size!

Versions

QGIS version 3.22.2-Białowieża QGIS code revision 1601ec46d0 Qt version 5.15.2 Python version 3.9.5 GDAL/OGR version 3.4.0 PROJ version 8.2.0 EPSG Registry database version v10.038 (2021-10-21) GEOS version 3.10.0-CAPI-1.16.0 SQLite version 3.35.2 PDAL version 2.3.0 PostgreSQL client version 13.0 SpatiaLite version 5.0.1 QWT version 6.1.3 QScintilla2 version 2.11.5 OS version Windows 10 Version 2009

Active Python plugins AnotherDXF2Shape 1.2.3 AutoLayoutTool 1.1 batchvectorlayersaver 0.9 BulkVectorExport 1.1 changeDataSource 3.1 DataPlotly 3.8.1 deactivate_active_labels 0.5 FlowEstimator 0.21 flowTrace 1.1.1 Generalizer3 1.0 GeoCoding 2.18 geo_sim_processing 1.2.0 getthemfiltered 0.1.3 gridSplitter 0.4.0 GroupStats 2.2.5 HideDocks 0.6.1 ImageServerConnector 2.1.1 ImportPhotos 3.0.3 joinmultiplelines Version 0.4.1 karika 1.5 LayerBoard 1.0.1 linz-data-importer 2.2.3 loadthemall 3.3.0 MagicWand-master 1.3.1 MemoryLayerSaver 4.0.4 mmqgis 2021.9.10 nominatim_locator_filter 0.2.4 numericalDigitize 0.4.6 pathfinder version 0.4.1 plaingeometryeditor 3.0.0 plugin_reloader 0.9.1 powerpan 2.0 processing_saga 0.5.0 processing_taudem 3.0.0 processing_whitebox 0.14.0 profiletool 4.2.1 qchainage 3.0.1 QCopycanvas 0.5 qgis-plugin-findreplace-main 1 Qgis2threejs 2.6 QGIS3-getWKT 1.4 QuickMultiAttributeEdit3 version 3.0.3 QuickPrint 3.6.1 quicksaveqml 0.1.5 rasmover-master version 0.2 redLayer 2.2 selectThemes 3.0.1 simple_tools 0.4.1 splitmultipart 1.0.0 SplitPolygonShowingAreas 0.13 statist 3.2 valuetool 3.0.15 ViewshedAnalysis 1.7 WaterNetAnalyzer-master 1.7 grassprovider 2.12.99 processing 2.12.99 sagaprovider 2.12.99 autoSaver 2.6 CalculateGeometry 0.6.4 numerator 0.2 quick_map_services 0.19.27 themeselector 3.2.2 volume_calculation_tool 0.4

Supported QGIS version

New profile

Additional context

No response

AlisterH commented 2 years ago

If I create a vrt file referring to the .asc file I've noticed some interesting behaviour which seems related:

C:\Processing>gdal_translate -a_srs EPSG:2105 test.asc test.vrt
Input file size is 95, 56

C:\Processing>type test.vrt
<VRTDataset rasterXSize="95" rasterYSize="56">
  <SRS dataAxisToSRSAxisMapping="2,1">PROJCS["NZGD2000 / Mount Eden 2000",GEOGCS["NZGD2000",DATUM["New_Zealand_Geodetic_Datum_2000",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","6167"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4167"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",-36.8797222222222],PARAMETER["central_meridian",174.764166666667],PARAMETER["scale_factor",0.9999],PARAMETER["false_easting",400000],PARAMETER["false_northing",800000],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Northing",NORTH],AXIS["Easting",EAST],AUTHORITY["EPSG","2105"]]</SRS>
  <GeoTransform>  4.1563000000000000e+05,  1.0000000000000000e+01,  0.0000000000000000e+00,  7.5938000000000000e+05,  0.0000000000000000e+00, -1.0000000000000000e+01</GeoTransform>
  <VRTRasterBand dataType="Float32" band="1" blockYSize="1">
    <Metadata>
      <MDI key="STATISTICS_MAXIMUM">78.616</MDI>
      <MDI key="STATISTICS_MEAN">71.720051459052</MDI>
      <MDI key="STATISTICS_MINIMUM">45.724</MDI>
      <MDI key="STATISTICS_STDDEV">5.7599521748039</MDI>
      <MDI key="STATISTICS_VALID_PERCENT">59.91</MDI>
    </Metadata>
    <NoDataValue>-9999</NoDataValue>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">test.asc</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="95" RasterYSize="56" DataType="Float32" BlockXSize="95" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="95" ySize="56" />
      <DstRect xOff="0" yOff="0" xSize="95" ySize="56" />
    </SimpleSource>
  </VRTRasterBand>
</VRTDataset>

Load the VRT file into QGIS and note that values are returned (e.g. in the valuetool plugin) to many decimal places. The VRT file is also edited when the layer is removed from QGIS (you may need to open the layer properties or something first); note the second instance of Float32 is changed to Float64:

C:\Processing>type test.vrt
<VRTDataset rasterXSize="95" rasterYSize="56">
  <SRS dataAxisToSRSAxisMapping="2,1">PROJCS["NZGD2000 / Mount Eden 2000",GEOGCS["NZGD2000",DATUM["New_Zealand_Geodetic_Datum_2000",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","6167"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4167"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",-36.8797222222222],PARAMETER["central_meridian",174.764166666667],PARAMETER["scale_factor",0.9999],PARAMETER["false_easting",400000],PARAMETER["false_northing",800000],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Northing",NORTH],AXIS["Easting",EAST],AUTHORITY["EPSG","2105"]]</SRS>
  <GeoTransform>  4.1563000000000000e+05,  1.0000000000000000e+01,  0.0000000000000000e+00,  7.5938000000000000e+05,  0.0000000000000000e+00, -1.0000000000000000e+01</GeoTransform>
  <VRTRasterBand dataType="Float32" band="1">
    <Metadata>
      <MDI key="STATISTICS_MAXIMUM">78.616</MDI>
      <MDI key="STATISTICS_MEAN">71.720051459052</MDI>
      <MDI key="STATISTICS_MINIMUM">45.724</MDI>
      <MDI key="STATISTICS_STDDEV">5.7599521748039</MDI>
      <MDI key="STATISTICS_VALID_PERCENT">59.91</MDI>
    </Metadata>
    <NoDataValue>-9999</NoDataValue>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">test.asc</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="95" RasterYSize="56" DataType="Float64" BlockXSize="95" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="95" ySize="56" />
      <DstRect xOff="0" yOff="0" xSize="95" ySize="56" />
    </SimpleSource>
  </VRTRasterBand>
</VRTDataset>

Loading the file into QGIS again, values are displayed to many decimal places still. But if I edit the VRT, changing the first instance of Float32 to Float64, they are displayed to 3dp (regardless of whether the second instance is Float32 or Float64). It makes sense to me for the first instance to be Float64, and the second instance Float32, but QGIS will always change the second instance to Float 64...