2

I am new to GIS and I want to merge xyz tiles to one raster tiff.

Here is the source to the xyz files: dgm1_05314000_Bonn_EPSG4647_XYZ.zip

I have downloaded the tiles as .xyz files, then sorted the coordinates with this bash command (since they are for some reason not in the standard format):

for file in *.xyz; 
    # sort coordinates
    do sort -k2 -n -k1 "$file" -o "$file";
done

gdalinfo for one of the tiles gives me:

Driver: XYZ/ASCII Gridded XYZ
Files: dgm1_32372_5622_2_nw.xyz
Size is 2000, 2000
Coordinate System is `'
Origin = (32371999.500000000000000,5621999.500000000000000)
Pixel Size = (1.000000000000000,1.000000000000000)
Corner Coordinates:
Upper Left  (32371999.500, 5621999.500) 
Lower Left  (32371999.500, 5623999.500) 
Upper Right (32373999.500, 5621999.500) 
Lower Right (32373999.500, 5623999.500) 
Center      (32372999.500, 5622999.500) 
Band 1 Block=2000x1 Type=Float32, ColorInterp=Undefined
  Min=87.800 Max=160.880 

Then I want to run to get a .vrt file:

gdalbuildvrt mosaic.vrt *.xyz

But I am getting the following error:

Warning 6: gdalbuildvrt does not support positive NS resolution. Skipping (filename).xyz

How can I now merge these tiles to one raster tiff?

Any alternatives are also welcome.

2 Answers 2

3

So here is a full outline of my solution. It seems that the xyz. tiles have an odd formatting since the coordinates are not in the expected order and have a positive N-S pixel resolution. To fix the coordinates I first sorted the coordinates:

for file in *.xyz; 
    # sort coordinates
    do sort -k2 -n -k1 "$file" -o "$file";
done

The positive N-S pixel resolution means that the origin of the images are not in the upper left corner of the image but in the lower left (which is quite unusual). To fix this issue I found a workaround here. The gdal warp command with the correct EPSG code fixes the issue:

for file in *.xyz; 
    # adjust NS resolution
    do echo "$file"; gdalwarp -t_srs EPSG:4647 -overwrite "$file" "$(basename "$file" .xyz).tif"
done

Then I could build a .vrt and convert it to a tif:

gdalbuildvrt mosaic.vrt *.tif
gdal_translate -of GTiff -a_srs EPSG:4647 mosaic.vrt  merged.tif
2
  • I think the sort is unnecessary, fwiw - you could also warp into VRT and save space but wth you solved it ;)
    – mdsumner
    Jan 22, 2020 at 15:44
  • I found gdal_translate to be slow and do the odd sorting. I posted my solution in Python, which is much faster, respects standard ordering and uses multiprocessing in this answer: gis.stackexchange.com/a/435515/208488
    – oekopez
    Jul 8, 2022 at 8:47
0

I had the same issues with this data. My approach was to create raster-tiles from the point data and then to mosaic it together.

However, I did it in R, and I am pretty sure there is a nicer solution, however it worked for me:

    library(devtools)
    if (!require(parallelsugar)) install_github('nathanvan/parallelsugar')
    library(parallelsugar)
    library(raster)
    library(sp)

    files <- list.files("dom1l_05515000_Münster_EPSG25832_XYZ",full.names = T,pattern = ".xyz$")

    mclapply(1:length(files), function(x){
      f <- files[x]
      pts <- read.table(f, header=FALSE, col.names=c("x", "y", "z"),sep = ",", dec = ".")
      coordinates(pts) = ~x+y
      rast <- raster(ext=extent(pts), resolution = 1)
      rasOut<-rasterize(pts, rast, pts$z, fun = mean)
      fout=paste0("dom1l_05515000_Münster_EPSG25832_XYZ/Tiles/Tiles_",x)
      writeRaster(rasOut, fout, format="GTiff")
      })

and then I built a mosaic out of the raster tiles

1
  • thanks for your quick answer but I found an alternative workaround
    – sehan2
    Jan 22, 2020 at 8:29

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.