gdal_translate creates images that are mirrored

Question

I'm having a problem with gdal_translate - when translating from an XYZ file to AAIGrid or geotiff the raster is mirrored and "moved" south on one axis.

Mirrored image, lidar data from some open mines

The top part of the attached image is the original file, when displayed in QGis, the bottom file has been translated to AAIGrid format with gdal_translate using

gdal_translate -of AAIGrid -a_srs EPSG:3301 grid.xyz grid.grd

The same happens when I try to translate the XYZ file to Geotiff.

Edit: added file sample. The file looks like this, the first 2 colums are coordinates in EPSG:3301 local coordinate system and the third column is height in metres

719935 6575005 30.709999
719945 6575005 31.08
719955 6575005 30.805
719965 6575005 30.772499
719975 6575005 30.2775
719985 6575005 30.1175
719995 6575005 30.012501
715005 6575015 28.525
715015 6575015 28.715
715025 6575015 28.834999
715035 6575015 28.6875
715045 6575015 28.452499
715055 6575015 29.147499

Output from gdalinfo

gdalinfo 65714_dem_10m.xyz.grid 
Driver: AAIGrid/Arc/Info ASCII Grid
Files: 65714_dem_10m.xyz.grid
       65714_dem_10m.xyz.prj
Size is 500, 500
Coordinate System is:
PROJCS["Estonian_Coordinate_System_of_1997",
    GEOGCS["GCS_EST97",
        DATUM["Estonia_1997",
            SPHEROID["GRS_1980",6378137,298.257222101]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.017453292519943295]],
    PROJECTION["Lambert_Conformal_Conic_2SP"],
    PARAMETER["standard_parallel_1",59.33333333333334],
    PARAMETER["standard_parallel_2",58],
    PARAMETER["latitude_of_origin",57.51755393055556],
    PARAMETER["central_meridian",24],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",6375000],
    UNIT["Meter",1]]
Origin = (715000.000000000000000,6575000.000000000000000)
Pixel Size = (10.000000000000000,-10.000000000000000)
Corner Coordinates:
Upper Left  (  715000.000, 6575000.000) ( 27d46'16.00"E, 59d15'32.02"N)
Lower Left  (  715000.000, 6570000.000) ( 27d45'58.29"E, 59d12'50.69"N)
Upper Right (  720000.000, 6575000.000) ( 27d51'31.03"E, 59d15'22.83"N)
Lower Right (  720000.000, 6570000.000) ( 27d51'12.91"E, 59d12'41.51"N)
Center      (  717500.000, 6572500.000) ( 27d48'44.56"E, 59d14' 6.79"N)
Band 1 Block=500x1 Type=Float32, ColorInterp=Undefined

It seems that XYZ file is written from south to north because Y is increasing 6575005 -> 6575015. That is OK for the GDAL driver http://www.gdal.org/frmt_xyz.html. However, I believe that flip is caused by that. What do you get with gdalinfo grid.xyz? — user30184, Feb 04 '15 at 14:33
It seems that GDAL does not understand correctly the write order of your XYZ file. It claims that the upper northing is 6575000 but in your data sample we can see bigger numbers like 6575015. Write a letter to gdal-dev mailing list about the issue. If your DEM is open data put one XYZ file available for download. — user30184, Feb 04 '15 at 15:06
Trying your sample data, I get the negative pixel size for the asc output, but the geotif made from the xyz (not the asc!) is ok. — AndreJ, Feb 04 '15 at 16:04

AndreJ · Accepted Answer · 2015-02-04T20:08:19.333

The .xyz file is sorted in an unusual way, X ascending and Y ascending. Most GDAL drivers can handle that orientation by flipping the image internally. Thus the "upper" and "lower" coordinates are flipped too:

Driver: XYZ/ASCII Gridded XYZ
Files: test.xyz
Size is 500, 2
Coordinate System is `'
Origin = (715000.000000000000000,6575000.000000000000000)
Pixel Size = (10.000000000000000,10.000000000000000)
Corner Coordinates:
Upper Left  (  715000.000, 6575000.000) 
Lower Left  (  715000.000, 6575020.000) 
Upper Right (  720000.000, 6575000.000) 
Lower Right (  720000.000, 6575020.000) 
Center      (  717500.000, 6575010.000) 
Band 1 Block=500x1 Type=Float32, ColorInterp=Undefined
  Min=28.452 Max=31.080 
  NoData Value=0

If you sort the .xyz file by X ascending and Y descending, the metadata appears in a logical way:

Driver: XYZ/ASCII Gridded XYZ
Files: testsort.xyz
Size is 500, 2
Coordinate System is `'
Origin = (715000.000000000000000,6575020.000000000000000)
Pixel Size = (10.000000000000000,-10.000000000000000)
Corner Coordinates:
Upper Left  (  715000.000, 6575020.000) 
Lower Left  (  715000.000, 6575000.000) 
Upper Right (  720000.000, 6575020.000) 
Lower Right (  720000.000, 6575000.000) 
Center      (  717500.000, 6575010.000) 
Band 1 Block=500x1 Type=Float32, ColorInterp=Undefined
  Min=28.452 Max=31.080 
  NoData Value=0

Unfortunately, the AAIGrid format is not designed for the first orientation, and noone has implemented to catch this issue by an error message or internally flipping back the coordinates. As a consequence, the header for AAIGrid gets generated wrong:

ncols        500
nrows        2
xllcorner    715000.000000000000
yllcorner    6574980.000000000000
cellsize     10.000000000000
NODATA_value  0

The following data is always linewise from upper left to lower right.

As a kind of hack, you can use gdalwarp instead of gdal_translate to save your raster as tif, and then translate to AAIGrid:

gdalwarp test.xyz 3301.tif
gdalinfo 3301.tif
gdal_translate -of AAIGRID 3301.tif 3301.asc

which delivers the correct header:

ncols        500
nrows        2
xllcorner    715000.000000000000
yllcorner    6575000.000000000000
cellsize     10.000000000000
NODATA_value  0

You are a life saver. I was hung up on this for hours. In my case I was generating a tif from gdal_grid then gdalwarp to a tif (with no arguments), finally calling translate to convert to AAIGRID. So annoying. — jocull, Mar 31 '16 at 15:59

score 1 · Answer 2 · answered Feb 04 '15 at 15:16

You can try if this workaround happens to work. Use this command:

gdal_translate -of GTiff -co PROFILE=BASELINE -co TFW=YES -a_srs EPSG:3301 grid.xyz grid.tif

Open the TFW file and edit the 4th line that gives the North-South pixel size. I believe that you have in now as "-10". Change it to positive "10" and test with QGIS. The image may open to a correct place. If it does, you can write a new GeoTIFF with all normal geotransformation tags with gdalwarp

gdalwarp -of GTiff -s_srs epsg:3301 -t_srs epsg:3301 grid.tif normalized.tif

I must say that this is just a wild guess but it may work.

score 0 · Answer 3 · answered May 05 '19 at 14:55

Resuming an old topic on a still-current issue.

I propose 2 workarounds to solve the vertical mirroring issue of ESRI ASC files exported from QGIS:

A GUI-only workaround
A command line workaround

GUI-only workaround

This entails:

Export as GeoTiff, while setting a negative "Vertical Resolution"
Use Raster->Conversion->Translate to export as ASC
Set ASC export options

Export as GeoTiff, while setting a negative "Vertical Resolution"

This creates a standard, valid GeoTIFF file, in which the order of rows is inverted. In order to do so, the "Vertical Resolution" needs to be negative. Just add a "minus" sign to whatever resolution is there. Do not change its absolute value. For example:

"Vertical Resolution = 1000" -> "Vertical Resolution = -1000"
"Vertical Resolution = 50" -> "Vertical Resolution = -50"

you get the idea :) Check the "Add saved file to map" in order to load the saved GeoTIFF file into your current map, and visually inspect it.

Use Raster->Conversion->Translate to export as ASC

ASC export in QGIS is only available trhough gdal_translate (and other GDAL tools, actually), which you can access via menu Raster->Conversion->Translate.

Select on the Layer panel the file which you just exported, and choose the menu command "Raster->Conversion->Translate":

Set ASC export options

Make sure to have these options checked. In particular, "Add saved file to map" will import the newly created asc into your current QGIS map so that you can visually make sure everything worked out correctly.

Click Run; if all was OK, your map should have the final ASC layer rendered exactly over the original raster.

Command line workaround

See the GIST I published: https://gist.github.com/rafdouglas/9a5548b4a259ec2c7823826b65feeafd

#!/bin/bash
#
# corrects the vertical mirror error of ASC files exported from QGIS
# through the "Raster->Conversion->Translate" command
# (which in turn invokes gdal_translate)
#
# RafDouglas - 190505
#
#
# example header section of the ASC file:
#
#ncols        600
#nrows        300
#xllcorner    565200.000000000000
#yllcorner    4797400.000000000000
#cellsize     10.000000000000
#NODATA_value  0

if [ -z $1 ]; then

    echo "Usage: $0 input_file.asc"
    exit
fi

if [ ! -e $1 ]; then

    echo "Usage: $0 input_file.asc"
    echo "The specified input file ($1) does not exist."
    exit
fi

mypid=$$

mytmp_file='temp_'"$mypid"'.asc'
mytmp_prj='temp_'"$mypid"'.prj'

infile=$1
inprj=$(echo $infile|sed 's/.asc/.prj/g')

outfile=$(echo $infile|sed 's/.asc/_proper.asc/g')
outprj=$(echo $infile|sed 's/.asc/_proper.prj/g')


head -6 "$infile"|awk '{h[$1]=$2}END{
yllcorner_proper=h["yllcorner"]+h["nrows"]*h["cellsize"]
print "ncols        " h["ncols"]
print "nrows        " h["nrows"]
print "xllcorner    " h["xllcorner"]
#print "yllcorner    " h["yllcorner"]
print "yllcorner    " yllcorner_proper
print "cellsize     " h["cellsize"]
print "NODATA_value " h["NODATA_value"]
}' > "$mytmp_file"

#add the data rows in reversed order 
tail -n+7 "$infile"|tac >> "$mytmp_file"

if [ ! -e $outfile ]; then
    mv -n "$mytmp_file" "$outfile"
    cp -n "$inprj" "$outprj"
    echo
    echo "Output file: $outfile is ready"
else
    cp -n "$inprj" "$mytmp_prj"
    echo
    echo "Waring: $outfile is already present. Not overwriting it."
    echo "Output file: $mytmp_file is ready"
fi

echo