PostGIS  3.0.6dev-r@@SVN_REVISION@@
ovdump.py
Go to the documentation of this file.
1 #! /usr/bin/env python
2 #
3 #
4 # This is a simple script based on GDAL to dump overview to separate file.
5 # It is used in WKTRaster testing to compare raster samples.
6 #
7 # NOTE: GDAL does include Frank's (unofficial) dumpoverviews utility too,
8 # which dumps overview as complete geospatial raster dataset.
9 #
10 # Copyright (C) 2009 Mateusz Loskot <mateusz@loskot.net>
11 #
12 # This program is free software; you can redistribute it and/or
13 # modify it under the terms of the GNU General Public License
14 # as published by the Free Software Foundation; either version 2
15 # of the License, or (at your option) any later version.
16 #
17 # This program is distributed in the hope that it will be useful,
18 # but WITHOUT ANY WARRANTY; without even the implied warranty of
19 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 # GNU General Public License for more details.
21 #
22 # You should have received a copy of the GNU General Public License
23 # along with this program; if not, write to the Free Software Foundation,
24 # Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
25 #
26 
27 from __future__ import print_function
28 from osgeo import gdal
29 from osgeo import osr
30 import osgeo.gdalconst as gdalc
31 from optparse import OptionParser
32 import os
33 import struct
34 import sys
35 
36 
37 prs = OptionParser(version="%prog $Revision: 4037 $",
38  usage="%prog -r <RASTER> -v <OVERVIEW>",
39  description="Dump GDAL raster overview to separate file, GeoTIF by default")
40 prs.add_option("-r", "--raster", dest="raster", action="store", default=None,
41  help="input raster file")
42 prs.add_option("-v", "--overview", dest="overview", action="store", type="int", default=None,
43  help="1-based index of raster overview, optional")
44 
45 (opts, args) = prs.parse_args()
46 
47 if opts.raster is None:
48  prs.error("use option -r to specify input raster file that contains overviews")
49 
50 
52 
53 # Source
54 src_ds = gdal.Open(opts.raster, gdalc.GA_ReadOnly);
55 if src_ds is None:
56  sys.exit('ERROR: Cannot open input file: ' + opts.raster)
57 
58 band = src_ds.GetRasterBand(1)
59 if opts.overview is None:
60  nvstart = 0
61  nvend = band.GetOverviewCount()
62 else:
63  nvstart = opts.overview - 1
64  nvend = opts.overview
65 
66 for nv in range(nvstart, nvend):
67 
68  band = src_ds.GetRasterBand(1)
69  if nv < 0 and nv >= band.GetOverviewCount():
70  print("ERROR: Failed to find overview %d" % nv)
71  sys.exit(1)
72  ov_band = band.GetOverview(nv)
73 
74  ovf = int(0.5 + band.XSize / float(ov_band.XSize))
75 
76  print("--- OVERVIEW #%d level = %d ---------------------------------------" % (nv + 1, ovf))
77 
78  # Destination
79  dst_file = os.path.basename(opts.raster) + "_ov_" + str(ovf) + ".tif"
80  dst_format = "GTiff"
81  dst_drv = gdal.GetDriverByName(dst_format)
82  dst_ds = dst_drv.Create(dst_file, ov_band.XSize, ov_band.YSize, src_ds.RasterCount, ov_band.DataType)
83 
84  print("Source: " + opts.raster)
85  print("Target: " + dst_file)
86  print("Exporting %d bands of size %d x %d" % (src_ds.RasterCount, ov_band.XSize, ov_band.YSize))
87 
88  # Rewrite bands of overview to output file
89  ov_band = None
90  band = None
91 
92  for nb in range(1, src_ds.RasterCount + 1):
93 
94  band = src_ds.GetRasterBand(nb)
95  assert band is not None
96  ov_band = band.GetOverview(nv)
97  assert ov_band is not None
98 
99  print(" Band #%d (%s) (%d x %d)" % \
100  (nb, gdal.GetDataTypeName(ov_band.DataType), ov_band.XSize, ov_band.YSize))
101 
102  dst_band = dst_ds.GetRasterBand(nb)
103  assert dst_band is not None
104  data = ov_band.ReadRaster(0, 0, ov_band.XSize, ov_band.YSize)
105  assert data is not None
106  dst_band.WriteRaster(0, 0, ov_band.XSize, ov_band.YSize, data, buf_type = ov_band.DataType)
107 
108  dst_ds = None
109 
110 # Cleanup
111 src_ds = None
#define str(s)