Hi, I do not know what happens, but it seems that the algorithm finds rather few points within the search ellipse. The radius of the circle is 40 m and by looking at the data there are typically at least 30 points within that distance. However, for filling the whole raster I had to use "min_points=4". I could repeat the behavior by opening the 3857 data into QGIS and running the invdist from there. And when I reprojected the point data into EPSG:32634 then QGIS created a good gridded raster. The difference in scale between EPSG:3857 and EPSG:32634 is not so big in Poland that it could be the reason IMHO. For removing possible errors due to CSV format and VRT I materialized the datasets into OpenJUMP JML format. It could have been any other GIS format.
Unfortunately I cannot help you more. If I were a programmer I think I would try to catch how many points the algorithm finds within the search radius in EPSG:3857, compare how many points it finds from EPSG:32634 data, and finally consider if the results make sense. -Jukka Rahkonen- Lähettäjä: Paul Meems <bontepaarden@gmail.c Lähetetty: keskiviikko 22. toukokuuta 2024 18.24 Vastaanottaja: gdal-dev@lists.osgeo.org Kopio: Rahkonen Jukka <jukka.rahko...@maanmittauslaitos.fi> Aihe: Re: [gdal-dev] Gdal_grid produces a lot of empty cells Sorry Jukka, Forgot all about that ;) I've attached the csv from the field in Poland with the original lon-lat values and the X and Y in EPSG:3857 Op wo 22 mei 2024 om 16:57 schreef Rahkonen Jukka <jukka.rahko...@maanmittauslaitos.fi<mailto:jukka.rahko...@maanmittauslaitos.fi>>: Hi, At least I cannot resolve your issue on paper. Please share some test data. -Jukka- Lähettäjä: gdal-dev <gdal-dev-boun...@lists.osgeo.org<mailto:gdal-dev-boun...@lists.osgeo.org>> Puolesta Paul Meems via gdal-dev Lähetetty: keskiviikko 22. toukokuuta 2024 17.31 Vastaanottaja: gdal-dev@lists.osgeo.org<mailto:gdal-dev@lists.osgeo.org> Aihe: [gdal-dev] Gdal_grid produces a lot of empty cells Hello List, I have a CSV-file with sensor data and GPS coordinates in WGS84. I convert these GPS coordinates to a projected projection. Mostly a UTM zone, depending on the location. These X and Y values are added to the CSV-file and saved with a different name. Then I create a VRT-file to get the proper Z-column: <OGRVRTDataSource> <OGRVRTLayer name="01 field data-EPSG32634"> <SrcDataSource>01 field data-EPSG32634.csv</SrcDataSource> <GeometryType>wkbPoint</GeometryType> <GeometryField encoding="PointFromColumns" x="X" y="Y" z="Countrate"/> </OGRVRTLayer> </OGRVRTDataSource> This VRT file is sent to Gdal_grid to create a tiff-file using these arguments: -a invdist:power=2: smoothing=20: min_points=16: max_points=56: radius1=40: radius2=40: nodata=1.701410000000000071e+38 -txe 273360 274940 -tye 5559440 5561020 -tr 12.5 12.5 -a_srs EPSG:32634 "01 field data-EPSG32634-TC.vrt" "tmpmrmq4d.tiff" This tiff-file looks like this: https://pasteboard.co/tKP1ZJ1cAY2r.png This process works fine and we've been using this for several years now. Recently we added a process step to export one of the sensor data to a PostGIS database and with the help of Geoserver and OpenLayers we've created a very simple webpage. In this step, things are getting weird. The WebMap is in EPSG:3857 (Pseudo/Web-Mercator). I start with the original CSV-file with the sensor data and add the X and Y columns again, now in EPSG:3857. Then I create a similar vrt-file: <OGRVRTDataSource> <OGRVRTLayer name="01 field data-EPSG3857"> <SrcDataSource>01 field data-EPSG3857.csv</SrcDataSource> <GeometryType>wkbPoint</GeometryType> <GeometryField encoding="PointFromColumns" x="X" y="Y" z="Countrate"/> </OGRVRTLayer> </OGRVRTDataSource> And call gdal_grid with similar arguments: -a invdist:power=2: smoothing=20: min_points=16: max_points=56: radius1=40: radius2=40: nodata=1.701410000000000071e+38 -txe 1984580 1987000 -tye 6471320 6473740 -tr 12.5 12.5 -a_srs EPSG:3857 "01 field data-EPSG3857-TC.vrt" "tmpun25fo.tiff" But now the result is: https://pasteboard.co/bcMp7rP5eyff.png A huge amount of empty cells. I don't understand this. Both EPSG:32634 (WGS 84 / UTM zone 34N) and EPSG:3857 are in meters so I thought the same settings should work. But in fact only with very different settings I can get a proper grid: -a invdist:power=2: smoothing=60: min_points=4: max_points=56: radius1=80: radius2=80: nodata=1.701410000000000071e+38 -txe 1984580 1987000 -tye 6471320 6473740 -tr 20 20 -a_srs EPSG:3857 "01 field data-EPSG3857-TC.vrt" "tmpsvlhgy.tiff" I got these values through trial and error. https://pasteboard.co/7pypGsynS09K.png But now we loose a lot of variation because the search radius and grid cell size are much higher. To make it even more weird. Other fields nearby don't have this problem and I also have this problem with other projections like: Belgium (EPSG:3857): https://pasteboard.co/ZXLTsEuvpZ3s.png Belgium (EPSG:31370): https://pasteboard.co/yNpF6yujF3tV.png The Netherlands (EPSG:3857): https://pasteboard.co/8PpDJrDi1Wxq.png The Netherlands (EPSG:28992): https://pasteboard.co/38SlnHVCog4m.png Norway (EPSG:3857): https://pasteboard.co/q69jQHRK40c6.png Norway (EPSG: 25833): https://pasteboard.co/gs9HNAmnaxSJ.png We do have enough data points: https://pasteboard.co/4CmBs7kgQGVO.png https://pasteboard.co/J8ZbMHKeErZu.png Thanks if you made it all the way here ;) My questions are: * Am I using Gdal_Grid in the proper way or do I need to change one or more arguments? * What might be an explanation for the empty cells when using EPSG:3857? Thanks, Paul Meems
_______________________________________________ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
