File: rtreader.py

package info (click to toggle)
postgis 3.5.2%2Bdfsg-1
  • links: PTS, VCS
  • area: main
  • in suites: sid, trixie
  • size: 70,052 kB
  • sloc: ansic: 162,204; sql: 93,950; xml: 53,121; cpp: 12,646; perl: 5,658; sh: 5,369; makefile: 3,434; python: 1,205; yacc: 447; lex: 151; pascal: 58
file content (191 lines) | stat: -rwxr-xr-x 7,016 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
#! /usr/bin/env python
#
#
# A simple driver to read RASTER field data directly from PostGIS/WKTRaster.
#
# Copyright (C) 2009 Mateusz Loskot <mateusz@loskot.net>
# 
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 2
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software Foundation,
# Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
#
###############################################################################
#
# Requirements: psycopg2
#
# WARNING: The main purpose of the RasterReader is to test and debug
# WKT Raster implementation. It is not supposed to be an efficient
# performance killer, by no means.
#
###############################################################################
import psycopg2
import sys

###############################################################################
# RASTER driver (read-only)

class RasterError(Exception):
    def __init__(self, msg):
        self.msg = msg 
    def __str__(self):
        return self.msg

class RasterReader(object):
    """Reader of RASTER data stored in specified column and row (where) in a table"""

    # Constructors

    def __init__(self, connstr, table, column, where = ""):
        self._connstr = connstr 
        self._conn = None
        self._table = table
        self._column = column
        self._where = where
        self._sizes = None
        self._types = None
        self._logging = False
        # Connect and read RASTER header
        self._setup()

    # Public properties

    logging = property(fset = lambda self, v: setattr(self, '_logging', v))
    db = property(fget = lambda self: self._get_db())
    table = property(fget = lambda self: self._table)
    column = property(fget = lambda self: self._column)
    width = property(fget = lambda self: self._get_width())
    height = property(fget = lambda self: self._get_height())
    num_bands = property(fget = lambda self: self._get_num_bands())
    pixel_types = property(fget = lambda self: self._get_pixel_types())

    # Public methods 

    def get_value(self, band, x, y):
        return self._query_value(band, x, y)

    def copy_to(self, file, raster_format='TIFF', output_format='HEX', sep='\t'):
        """
        Proxy for SQL command COPY TO,
        Converts selected rasters to specified raster_format with output sent either to 
        single hex-based plain text file or one or more binary files in raster_format,
        one raster binary file per tuple from the raster table.
        The BIN output uses HEX output as intermediate stage.
        raster_format - TIFF|JPEG|PNG
        output_format - HEX|BIN; BIN is a binary file in raster_format
        sep - if output_format=HEX, separates rid value from hex-encoded binary.
        """
        import os.path
        filehex = file # No extension added, may be user-defined
        with open(filehex, 'w') as f:
            select = "SELECT rid, encode(ST_As%s(%s), 'hex') As rt FROM %s" % (raster_format, self._column, self._table)
            if self._where is not None and len(self._where) > 0:
                select += ' WHERE %s' % self._where
            sql = "COPY (%s) TO STDOUT (DELIMITER '%s')" % (select, sep)        
            cur = self._conn.cursor()
            cur.copy_expert(sql, f)

        if output_format == 'BIN':
            import binascii
            with open(filehex, 'r') as f:
                dirname = os.path.dirname(file)
                ext = raster_format.lower()
                for line in f.readlines():
                    rid, raster = line.split()
                    filebin = self._table + '_' + self._column + '_' + rid + '.' + ext
                    filebin = os.path.join(dirname, filebin)
                    with open(filebin, 'w+') as fbin:
                        fbin.write(binascii.unhexlify(raster))

    # Private methods

    def _log(self, m):
        if self._logging:
            sys.stderr.write('[rtreader] ' + str(m) + '\n')

    def _get_db(self):
        n = filter(lambda db: db[:6] == 'dbname', self._connstr.split())[0].split('=')[1]
        return n.strip('\'').strip()

    def _get_width(self):
        return self._query_raster_size(0)

    def _get_height(self):
        return self._query_raster_size(1)

    def _get_num_bands(self):
        return self._query_raster_size(2)

    def _get_pixel_types(self):
        return self._query_pixel_types()

    def _setup(self):
        self._connect()

    def _connect(self):
        try:
            if self._conn is None:
                self._conn = psycopg2.connect(self._connstr)
        except Exception as e:
            raise RasterError("Failed to connect to %s: %s" % (self._connstr, e))

    def _query_single_row(self, sql):
        assert self._conn is not None
        #self._log(sql)

        try:
            cur = self._conn.cursor()
            cur.execute(sql)
        except Exception as e:
            raise RasterError("Failed to execute query %s: %s" % (sql, e))

        row = cur.fetchone()
        if row is None:
            raise RasterError("No tupes returned for query: %s" % sql)
        return row

    def _query_value(self, band, x, y):
        sql = 'SELECT st_value(%s, %d, %d, %d) FROM %s' % \
                 (self._column, band, x, y, self._table)
        if len(self._where) > 0:
            sql += ' WHERE %s' % self._where

        row = self._query_single_row(sql)
        if row is None:
            raise RasterError("Value of pixel %dx%d of band %d is none" %(x, y, band))
        return row[0]
    
    def _query_raster_size(self, dim, force = False):
        if self._sizes is None or force is True:
            sql = 'SELECT st_width(%s), st_height(%s), st_numbands(%s) FROM %s' % \
                     (self._column, self._column, self._column, self._table)
            if len(self._where) > 0:
                sql += ' WHERE %s' % self._where
                
            self._log(sql)
            self._sizes = self._query_single_row(sql)

        if self._sizes is None:
            raise RasterError("Failed to query raster size of dim {} with force {}".format(dim, force))
        return self._sizes[dim]

    def _query_pixel_types(self):

        types = []
        sql = 'SELECT '
        for i in range(0, self.num_bands):
            if i != 0:
                sql += ','
            nband = i + 1
            sql += ' st_bandpixeltype(%s, %d) ' % (self._column, nband)
        sql += ' FROM ' + self._table
        return self._query_single_row(sql)