Salut, Almost there, it looks like the data is not inverted, it still starts from the upper left corner instead of lower left.
Thanks Denis On Sat, Nov 9, 2013 at 8:46 PM, Even Rouault <[email protected]>wrote: > Le samedi 09 novembre 2013 20:39:50, denis cohen a écrit : > > Hi, > > > > Thanks for the script. It works fine except for the actual data > (elevation > > in my case, this is a DEM). > > It does not write the correct values. > > Ah right, this version only worked for Byte data type. The below version > should work for any GDAL data type > > ----- Begin of script invert.py ------ > > # -*- coding: utf-8 -*- > import sys > from osgeo import gdal > > src_ds = gdal.Open(sys.argv[1]) > xsize = src_ds.RasterXSize > ysize = src_ds.RasterYSize > bandcount = src_ds.RasterCount > datatype = src_ds.GetRasterBand(1).DataType > out_drv = gdal.GetDriverByName('GTiff') > out_ds = out_drv.Create(sys.argv[2], xsize, ysize, bandcount, datatype) > out_ds.SetProjection(src_ds.GetProjectionRef()) > src_gt = src_ds.GetGeoTransform() > out_gt = [ src_gt[i] for i in range(6) ] > out_gt[3] = src_gt[3] + ysize * src_gt[5] > out_gt[5] = -src_gt[5] > out_ds.SetGeoTransform(out_gt) > for j in range(src_ds.RasterYSize): > data = src_ds.ReadRaster(0, ysize - 1 - j, xsize, 1, buf_type = > datatype) > out_ds.WriteRaster(0, j, xsize, 1, data, buf_type = datatype) > out_ds = None > src_ds = None > > ----- End of script ------ > > > > > > > Merci > > > > Denis > > > > > > > > On Sat, Nov 9, 2013 at 8:02 PM, Even Rouault > > > > <[email protected]>wrote: > > > Le samedi 09 novembre 2013 10:46:20, denis cohen a écrit : > > > > Hello, > > > > > > > > I am a new user of gdal. > > > > I used gdal_translate to convert a netcdf file to xyz: > > > > gdal_translate -of XYZ file.nc file.xyz > > > > > > > > However, the xyz file starts with the upperleft coordinate. > > > > Is it possible to have it write the xyz file starting with the lower > > > > left coordinate? > > > > > > Not directly although it wouldn't be complicated to add an option to do > > > that. > > > > > > You can try the following GDAL Python script that will create an > > > intermediary > > > TIFF image that has both inverted georeferencing and imagery in a > > > consistant > > > way. > > > > > > python invert.py file.nc file_reserved.tif > > > gdal_translate -of XYZ file_reversed.tif file.xyz > > > > > > > > > ----- Begin of script invert.py ------ > > > > > > # -*- coding: utf-8 -*- > > > import sys > > > from osgeo import gdal > > > > > > src_ds = gdal.Open(sys.argv[1]) > > > xsize = src_ds.RasterXSize > > > ysize = src_ds.RasterYSize > > > bandcount = src_ds.RasterCount > > > out_drv = gdal.GetDriverByName('GTiff') > > > out_ds = out_drv.Create(sys.argv[2], xsize, ysize, bandcount) > > > out_ds.SetProjection(src_ds.GetProjectionRef()) > > > src_gt = src_ds.GetGeoTransform() > > > out_gt = [ src_gt[i] for i in range(6) ] > > > out_gt[3] = src_gt[3] + ysize * src_gt[5] > > > out_gt[5] = -src_gt[5] > > > out_ds.SetGeoTransform(out_gt) > > > > > > for j in range(src_ds.RasterYSize): > > > data = src_ds.ReadRaster(0, ysize - 1 - j, xsize, 1) > > > out_ds.WriteRaster(0, j, xsize, 1, data) > > > > > > out_ds = None > > > src_ds = None > > > > > > ----- End of script ------ > > > > > > Even > > > > > > -- > > > Geospatial professional services > > > http://even.rouault.free.fr/services.html > > -- > Geospatial professional services > http://even.rouault.free.fr/services.html >
_______________________________________________ gdal-dev mailing list [email protected] http://lists.osgeo.org/mailman/listinfo/gdal-dev
