Ajoutez des fichiers projet.

This commit is contained in:
Ambulance Clerc
2021-12-18 18:43:17 +01:00
parent 3c4d48ed26
commit 46254605fc
4842 changed files with 732322 additions and 0 deletions

View File

@@ -0,0 +1,69 @@
"""
This object provides quoting for GEOS geometries into PostgreSQL/PostGIS.
"""
from psycopg2 import Binary
from psycopg2.extensions import ISQLQuote
from django.contrib.gis.db.backends.postgis.pgraster import to_pgraster
from django.contrib.gis.geos import GEOSGeometry
class PostGISAdapter:
def __init__(self, obj, geography=False):
"""
Initialize on the spatial object.
"""
self.is_geometry = isinstance(obj, (GEOSGeometry, PostGISAdapter))
# Getting the WKB (in string form, to allow easy pickling of
# the adaptor) and the SRID from the geometry or raster.
if self.is_geometry:
self.ewkb = bytes(obj.ewkb)
self._adapter = Binary(self.ewkb)
else:
self.ewkb = to_pgraster(obj)
self.srid = obj.srid
self.geography = geography
def __conform__(self, proto):
"""Does the given protocol conform to what Psycopg2 expects?"""
if proto == ISQLQuote:
return self
else:
raise Exception('Error implementing psycopg2 protocol. Is psycopg2 installed?')
def __eq__(self, other):
return isinstance(other, PostGISAdapter) and self.ewkb == other.ewkb
def __hash__(self):
return hash(self.ewkb)
def __str__(self):
return self.getquoted()
@classmethod
def _fix_polygon(cls, poly):
return poly
def prepare(self, conn):
"""
This method allows escaping the binary in the style required by the
server's `standard_conforming_string` setting.
"""
if self.is_geometry:
self._adapter.prepare(conn)
def getquoted(self):
"""
Return a properly quoted string for use in PostgreSQL/PostGIS.
"""
if self.is_geometry:
# Psycopg will figure out whether to use E'\\000' or '\000'.
return '%s(%s)' % (
'ST_GeogFromWKB' if self.geography else 'ST_GeomFromEWKB',
self._adapter.getquoted().decode()
)
else:
# For rasters, add explicit type cast to WKB string.
return "'%s'::raster" % self.ewkb

View File

@@ -0,0 +1,26 @@
from django.db.backends.base.base import NO_DB_ALIAS
from django.db.backends.postgresql.base import (
DatabaseWrapper as Psycopg2DatabaseWrapper,
)
from .features import DatabaseFeatures
from .introspection import PostGISIntrospection
from .operations import PostGISOperations
from .schema import PostGISSchemaEditor
class DatabaseWrapper(Psycopg2DatabaseWrapper):
SchemaEditorClass = PostGISSchemaEditor
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
if kwargs.get('alias', '') != NO_DB_ALIAS:
self.features = DatabaseFeatures(self)
self.ops = PostGISOperations(self)
self.introspection = PostGISIntrospection(self)
def prepare_database(self):
super().prepare_database()
# Check that postgis extension is installed.
with self.cursor() as cursor:
cursor.execute("CREATE EXTENSION IF NOT EXISTS postgis")

View File

@@ -0,0 +1,52 @@
"""
PostGIS to GDAL conversion constant definitions
"""
# Lookup to convert pixel type values from GDAL to PostGIS
GDAL_TO_POSTGIS = [None, 4, 6, 5, 8, 7, 10, 11, None, None, None, None]
# Lookup to convert pixel type values from PostGIS to GDAL
POSTGIS_TO_GDAL = [1, 1, 1, 3, 1, 3, 2, 5, 4, None, 6, 7, None, None]
# Struct pack structure for raster header, the raster header has the
# following structure:
#
# Endianness, PostGIS raster version, number of bands, scale, origin,
# skew, srid, width, and height.
#
# Scale, origin, and skew have x and y values. PostGIS currently uses
# a fixed endianness (1) and there is only one version (0).
POSTGIS_HEADER_STRUCTURE = 'B H H d d d d d d i H H'
# Lookup values to convert GDAL pixel types to struct characters. This is
# used to pack and unpack the pixel values of PostGIS raster bands.
GDAL_TO_STRUCT = [
None, 'B', 'H', 'h', 'L', 'l', 'f', 'd',
None, None, None, None,
]
# Size of the packed value in bytes for different numerical types.
# This is needed to cut chunks of band data out of PostGIS raster strings
# when decomposing them into GDALRasters.
# See https://docs.python.org/library/struct.html#format-characters
STRUCT_SIZE = {
'b': 1, # Signed char
'B': 1, # Unsigned char
'?': 1, # _Bool
'h': 2, # Short
'H': 2, # Unsigned short
'i': 4, # Integer
'I': 4, # Unsigned Integer
'l': 4, # Long
'L': 4, # Unsigned Long
'f': 4, # Float
'd': 8, # Double
}
# Pixel type specifies type of pixel values in a band. Storage flag specifies
# whether the band data is stored as part of the datum or is to be found on the
# server's filesystem. There are currently 11 supported pixel value types, so 4
# bits are enough to account for all. Reserve the upper 4 bits for generic
# flags.
# See https://trac.osgeo.org/postgis/wiki/WKTRaster/RFC/RFC1_V0SerialFormat#Pixeltypeandstorageflag
BANDTYPE_PIXTYPE_MASK = 0x0F
BANDTYPE_FLAG_HASNODATA = 1 << 6

View File

@@ -0,0 +1,13 @@
from django.contrib.gis.db.backends.base.features import BaseSpatialFeatures
from django.db.backends.postgresql.features import (
DatabaseFeatures as Psycopg2DatabaseFeatures,
)
class DatabaseFeatures(BaseSpatialFeatures, Psycopg2DatabaseFeatures):
supports_geography = True
supports_3d_storage = True
supports_3d_functions = True
supports_raster = True
supports_empty_geometries = True
empty_intersection_returns_none = False

View File

@@ -0,0 +1,60 @@
from django.contrib.gis.gdal import OGRGeomType
from django.db.backends.postgresql.introspection import DatabaseIntrospection
class PostGISIntrospection(DatabaseIntrospection):
postgis_oid_lookup = {} # Populated when introspection is performed.
ignored_tables = DatabaseIntrospection.ignored_tables + [
'geography_columns',
'geometry_columns',
'raster_columns',
'spatial_ref_sys',
'raster_overviews',
]
def get_field_type(self, data_type, description):
if not self.postgis_oid_lookup:
# Query PostgreSQL's pg_type table to determine the OID integers
# for the PostGIS data types used in reverse lookup (the integers
# may be different across versions). To prevent unnecessary
# requests upon connection initialization, the `data_types_reverse`
# dictionary isn't updated until introspection is performed here.
with self.connection.cursor() as cursor:
cursor.execute("SELECT oid, typname FROM pg_type WHERE typname IN ('geometry', 'geography')")
self.postgis_oid_lookup = dict(cursor.fetchall())
self.data_types_reverse.update((oid, 'GeometryField') for oid in self.postgis_oid_lookup)
return super().get_field_type(data_type, description)
def get_geometry_type(self, table_name, description):
"""
The geometry type OID used by PostGIS does not indicate the particular
type of field that a geometry column is (e.g., whether it's a
PointField or a PolygonField). Thus, this routine queries the PostGIS
metadata tables to determine the geometry type.
"""
with self.connection.cursor() as cursor:
cursor.execute("""
SELECT t.coord_dimension, t.srid, t.type FROM (
SELECT * FROM geometry_columns
UNION ALL
SELECT * FROM geography_columns
) AS t WHERE t.f_table_name = %s AND t.f_geometry_column = %s
""", (table_name, description.name))
row = cursor.fetchone()
if not row:
raise Exception('Could not find a geometry or geography column for "%s"."%s"' %
(table_name, description.name))
dim, srid, field_type = row
# OGRGeomType does not require GDAL and makes it easy to convert
# from OGC geom type name to Django field.
field_type = OGRGeomType(field_type).django
# Getting any GeometryField keyword arguments that are not the default.
field_params = {}
if self.postgis_oid_lookup.get(description.type_code) == 'geography':
field_params['geography'] = True
if srid != 4326:
field_params['srid'] = srid
if dim != 2:
field_params['dim'] = dim
return field_type, field_params

View File

@@ -0,0 +1,70 @@
"""
The GeometryColumns and SpatialRefSys models for the PostGIS backend.
"""
from django.contrib.gis.db.backends.base.models import SpatialRefSysMixin
from django.db import models
class PostGISGeometryColumns(models.Model):
"""
The 'geometry_columns' view from PostGIS. See the PostGIS
documentation at Ch. 4.3.2.
"""
f_table_catalog = models.CharField(max_length=256)
f_table_schema = models.CharField(max_length=256)
f_table_name = models.CharField(max_length=256)
f_geometry_column = models.CharField(max_length=256)
coord_dimension = models.IntegerField()
srid = models.IntegerField(primary_key=True)
type = models.CharField(max_length=30)
class Meta:
app_label = 'gis'
db_table = 'geometry_columns'
managed = False
def __str__(self):
return '%s.%s - %dD %s field (SRID: %d)' % (
self.f_table_name,
self.f_geometry_column,
self.coord_dimension,
self.type,
self.srid,
)
@classmethod
def table_name_col(cls):
"""
Return the name of the metadata column used to store the feature table
name.
"""
return 'f_table_name'
@classmethod
def geom_col_name(cls):
"""
Return the name of the metadata column used to store the feature
geometry column.
"""
return 'f_geometry_column'
class PostGISSpatialRefSys(models.Model, SpatialRefSysMixin):
"""
The 'spatial_ref_sys' table from PostGIS. See the PostGIS
documentation at Ch. 4.2.1.
"""
srid = models.IntegerField(primary_key=True)
auth_name = models.CharField(max_length=256)
auth_srid = models.IntegerField()
srtext = models.CharField(max_length=2048)
proj4text = models.CharField(max_length=2048)
class Meta:
app_label = 'gis'
db_table = 'spatial_ref_sys'
managed = False
@property
def wkt(self):
return self.srtext

View File

@@ -0,0 +1,390 @@
import re
from django.conf import settings
from django.contrib.gis.db.backends.base.operations import (
BaseSpatialOperations,
)
from django.contrib.gis.db.backends.utils import SpatialOperator
from django.contrib.gis.db.models import GeometryField, RasterField
from django.contrib.gis.gdal import GDALRaster
from django.contrib.gis.geos.geometry import GEOSGeometryBase
from django.contrib.gis.geos.prototypes.io import wkb_r
from django.contrib.gis.measure import Distance
from django.core.exceptions import ImproperlyConfigured
from django.db import NotSupportedError, ProgrammingError
from django.db.backends.postgresql.operations import DatabaseOperations
from django.db.models import Func, Value
from django.utils.functional import cached_property
from django.utils.version import get_version_tuple
from .adapter import PostGISAdapter
from .models import PostGISGeometryColumns, PostGISSpatialRefSys
from .pgraster import from_pgraster
# Identifier to mark raster lookups as bilateral.
BILATERAL = 'bilateral'
class PostGISOperator(SpatialOperator):
def __init__(self, geography=False, raster=False, **kwargs):
# Only a subset of the operators and functions are available for the
# geography type.
self.geography = geography
# Only a subset of the operators and functions are available for the
# raster type. Lookups that don't support raster will be converted to
# polygons. If the raster argument is set to BILATERAL, then the
# operator cannot handle mixed geom-raster lookups.
self.raster = raster
super().__init__(**kwargs)
def as_sql(self, connection, lookup, template_params, *args):
if lookup.lhs.output_field.geography and not self.geography:
raise ValueError('PostGIS geography does not support the "%s" '
'function/operator.' % (self.func or self.op,))
template_params = self.check_raster(lookup, template_params)
return super().as_sql(connection, lookup, template_params, *args)
def check_raster(self, lookup, template_params):
spheroid = lookup.rhs_params and lookup.rhs_params[-1] == 'spheroid'
# Check which input is a raster.
lhs_is_raster = lookup.lhs.field.geom_type == 'RASTER'
rhs_is_raster = isinstance(lookup.rhs, GDALRaster)
# Look for band indices and inject them if provided.
if lookup.band_lhs is not None and lhs_is_raster:
if not self.func:
raise ValueError('Band indices are not allowed for this operator, it works on bbox only.')
template_params['lhs'] = '%s, %s' % (template_params['lhs'], lookup.band_lhs)
if lookup.band_rhs is not None and rhs_is_raster:
if not self.func:
raise ValueError('Band indices are not allowed for this operator, it works on bbox only.')
template_params['rhs'] = '%s, %s' % (template_params['rhs'], lookup.band_rhs)
# Convert rasters to polygons if necessary.
if not self.raster or spheroid:
# Operators without raster support.
if lhs_is_raster:
template_params['lhs'] = 'ST_Polygon(%s)' % template_params['lhs']
if rhs_is_raster:
template_params['rhs'] = 'ST_Polygon(%s)' % template_params['rhs']
elif self.raster == BILATERAL:
# Operators with raster support but don't support mixed (rast-geom)
# lookups.
if lhs_is_raster and not rhs_is_raster:
template_params['lhs'] = 'ST_Polygon(%s)' % template_params['lhs']
elif rhs_is_raster and not lhs_is_raster:
template_params['rhs'] = 'ST_Polygon(%s)' % template_params['rhs']
return template_params
class ST_Polygon(Func):
function = 'ST_Polygon'
def __init__(self, expr):
super().__init__(expr)
expr = self.source_expressions[0]
if isinstance(expr, Value) and not expr._output_field_or_none:
self.source_expressions[0] = Value(expr.value, output_field=RasterField(srid=expr.value.srid))
@cached_property
def output_field(self):
return GeometryField(srid=self.source_expressions[0].field.srid)
class PostGISOperations(BaseSpatialOperations, DatabaseOperations):
name = 'postgis'
postgis = True
geom_func_prefix = 'ST_'
Adapter = PostGISAdapter
collect = geom_func_prefix + 'Collect'
extent = geom_func_prefix + 'Extent'
extent3d = geom_func_prefix + '3DExtent'
length3d = geom_func_prefix + '3DLength'
makeline = geom_func_prefix + 'MakeLine'
perimeter3d = geom_func_prefix + '3DPerimeter'
unionagg = geom_func_prefix + 'Union'
gis_operators = {
'bbcontains': PostGISOperator(op='~', raster=True),
'bboverlaps': PostGISOperator(op='&&', geography=True, raster=True),
'contained': PostGISOperator(op='@', raster=True),
'overlaps_left': PostGISOperator(op='&<', raster=BILATERAL),
'overlaps_right': PostGISOperator(op='&>', raster=BILATERAL),
'overlaps_below': PostGISOperator(op='&<|'),
'overlaps_above': PostGISOperator(op='|&>'),
'left': PostGISOperator(op='<<'),
'right': PostGISOperator(op='>>'),
'strictly_below': PostGISOperator(op='<<|'),
'strictly_above': PostGISOperator(op='|>>'),
'same_as': PostGISOperator(op='~=', raster=BILATERAL),
'exact': PostGISOperator(op='~=', raster=BILATERAL), # alias of same_as
'contains': PostGISOperator(func='ST_Contains', raster=BILATERAL),
'contains_properly': PostGISOperator(func='ST_ContainsProperly', raster=BILATERAL),
'coveredby': PostGISOperator(func='ST_CoveredBy', geography=True, raster=BILATERAL),
'covers': PostGISOperator(func='ST_Covers', geography=True, raster=BILATERAL),
'crosses': PostGISOperator(func='ST_Crosses'),
'disjoint': PostGISOperator(func='ST_Disjoint', raster=BILATERAL),
'equals': PostGISOperator(func='ST_Equals'),
'intersects': PostGISOperator(func='ST_Intersects', geography=True, raster=BILATERAL),
'overlaps': PostGISOperator(func='ST_Overlaps', raster=BILATERAL),
'relate': PostGISOperator(func='ST_Relate'),
'touches': PostGISOperator(func='ST_Touches', raster=BILATERAL),
'within': PostGISOperator(func='ST_Within', raster=BILATERAL),
'dwithin': PostGISOperator(func='ST_DWithin', geography=True, raster=BILATERAL),
}
unsupported_functions = set()
select = '%s::bytea'
select_extent = None
@cached_property
def function_names(self):
function_names = {
'AsWKB': 'ST_AsBinary',
'AsWKT': 'ST_AsText',
'BoundingCircle': 'ST_MinimumBoundingCircle',
'NumPoints': 'ST_NPoints',
}
if self.spatial_version < (2, 4, 0):
function_names['ForcePolygonCW'] = 'ST_ForceRHR'
return function_names
@cached_property
def spatial_version(self):
"""Determine the version of the PostGIS library."""
# Trying to get the PostGIS version because the function
# signatures will depend on the version used. The cost
# here is a database query to determine the version, which
# can be mitigated by setting `POSTGIS_VERSION` with a 3-tuple
# comprising user-supplied values for the major, minor, and
# subminor revision of PostGIS.
if hasattr(settings, 'POSTGIS_VERSION'):
version = settings.POSTGIS_VERSION
else:
# Run a basic query to check the status of the connection so we're
# sure we only raise the error below if the problem comes from
# PostGIS and not from PostgreSQL itself (see #24862).
self._get_postgis_func('version')
try:
vtup = self.postgis_version_tuple()
except ProgrammingError:
raise ImproperlyConfigured(
'Cannot determine PostGIS version for database "%s" '
'using command "SELECT postgis_lib_version()". '
'GeoDjango requires at least PostGIS version 2.4. '
'Was the database created from a spatial database '
'template?' % self.connection.settings_dict['NAME']
)
version = vtup[1:]
return version
def convert_extent(self, box):
"""
Return a 4-tuple extent for the `Extent` aggregate by converting
the bounding box text returned by PostGIS (`box` argument), for
example: "BOX(-90.0 30.0, -85.0 40.0)".
"""
if box is None:
return None
ll, ur = box[4:-1].split(',')
xmin, ymin = map(float, ll.split())
xmax, ymax = map(float, ur.split())
return (xmin, ymin, xmax, ymax)
def convert_extent3d(self, box3d):
"""
Return a 6-tuple extent for the `Extent3D` aggregate by converting
the 3d bounding-box text returned by PostGIS (`box3d` argument), for
example: "BOX3D(-90.0 30.0 1, -85.0 40.0 2)".
"""
if box3d is None:
return None
ll, ur = box3d[6:-1].split(',')
xmin, ymin, zmin = map(float, ll.split())
xmax, ymax, zmax = map(float, ur.split())
return (xmin, ymin, zmin, xmax, ymax, zmax)
def geo_db_type(self, f):
"""
Return the database field type for the given spatial field.
"""
if f.geom_type == 'RASTER':
return 'raster'
# Type-based geometries.
# TODO: Support 'M' extension.
if f.dim == 3:
geom_type = f.geom_type + 'Z'
else:
geom_type = f.geom_type
if f.geography:
if f.srid != 4326:
raise NotSupportedError('PostGIS only supports geography columns with an SRID of 4326.')
return 'geography(%s,%d)' % (geom_type, f.srid)
else:
return 'geometry(%s,%d)' % (geom_type, f.srid)
def get_distance(self, f, dist_val, lookup_type):
"""
Retrieve the distance parameters for the given geometry field,
distance lookup value, and the distance lookup type.
This is the most complex implementation of the spatial backends due to
what is supported on geodetic geometry columns vs. what's available on
projected geometry columns. In addition, it has to take into account
the geography column type.
"""
# Getting the distance parameter
value = dist_val[0]
# Shorthand boolean flags.
geodetic = f.geodetic(self.connection)
geography = f.geography
if isinstance(value, Distance):
if geography:
dist_param = value.m
elif geodetic:
if lookup_type == 'dwithin':
raise ValueError('Only numeric values of degree units are '
'allowed on geographic DWithin queries.')
dist_param = value.m
else:
dist_param = getattr(value, Distance.unit_attname(f.units_name(self.connection)))
else:
# Assuming the distance is in the units of the field.
dist_param = value
return [dist_param]
def get_geom_placeholder(self, f, value, compiler):
"""
Provide a proper substitution value for Geometries or rasters that are
not in the SRID of the field. Specifically, this routine will
substitute in the ST_Transform() function call.
"""
transform_func = self.spatial_function_name('Transform')
if hasattr(value, 'as_sql'):
if value.field.srid == f.srid:
placeholder = '%s'
else:
placeholder = '%s(%%s, %s)' % (transform_func, f.srid)
return placeholder
# Get the srid for this object
if value is None:
value_srid = None
else:
value_srid = value.srid
# Adding Transform() to the SQL placeholder if the value srid
# is not equal to the field srid.
if value_srid is None or value_srid == f.srid:
placeholder = '%s'
else:
placeholder = '%s(%%s, %s)' % (transform_func, f.srid)
return placeholder
def _get_postgis_func(self, func):
"""
Helper routine for calling PostGIS functions and returning their result.
"""
# Close out the connection. See #9437.
with self.connection.temporary_connection() as cursor:
cursor.execute('SELECT %s()' % func)
return cursor.fetchone()[0]
def postgis_geos_version(self):
"Return the version of the GEOS library used with PostGIS."
return self._get_postgis_func('postgis_geos_version')
def postgis_lib_version(self):
"Return the version number of the PostGIS library used with PostgreSQL."
return self._get_postgis_func('postgis_lib_version')
def postgis_proj_version(self):
"""Return the version of the PROJ library used with PostGIS."""
return self._get_postgis_func('postgis_proj_version')
def postgis_version(self):
"Return PostGIS version number and compile-time options."
return self._get_postgis_func('postgis_version')
def postgis_full_version(self):
"Return PostGIS version number and compile-time options."
return self._get_postgis_func('postgis_full_version')
def postgis_version_tuple(self):
"""
Return the PostGIS version as a tuple (version string, major,
minor, subminor).
"""
version = self.postgis_lib_version()
return (version,) + get_version_tuple(version)
def proj_version_tuple(self):
"""
Return the version of PROJ used by PostGIS as a tuple of the
major, minor, and subminor release numbers.
"""
proj_regex = re.compile(r'(\d+)\.(\d+)\.(\d+)')
proj_ver_str = self.postgis_proj_version()
m = proj_regex.search(proj_ver_str)
if m:
return tuple(map(int, m.groups()))
else:
raise Exception('Could not determine PROJ version from PostGIS.')
def spatial_aggregate_name(self, agg_name):
if agg_name == 'Extent3D':
return self.extent3d
else:
return self.geom_func_prefix + agg_name
# Routines for getting the OGC-compliant models.
def geometry_columns(self):
return PostGISGeometryColumns
def spatial_ref_sys(self):
return PostGISSpatialRefSys
def parse_raster(self, value):
"""Convert a PostGIS HEX String into a dict readable by GDALRaster."""
return from_pgraster(value)
def distance_expr_for_lookup(self, lhs, rhs, **kwargs):
return super().distance_expr_for_lookup(
self._normalize_distance_lookup_arg(lhs),
self._normalize_distance_lookup_arg(rhs),
**kwargs
)
@staticmethod
def _normalize_distance_lookup_arg(arg):
is_raster = (
arg.field.geom_type == 'RASTER'
if hasattr(arg, 'field') else
isinstance(arg, GDALRaster)
)
return ST_Polygon(arg) if is_raster else arg
def get_geometry_converter(self, expression):
read = wkb_r().read
geom_class = expression.output_field.geom_class
def converter(value, expression, connection):
return None if value is None else GEOSGeometryBase(read(value), geom_class)
return converter
def get_area_att_for_field(self, field):
return 'sq_m'

View File

@@ -0,0 +1,138 @@
import struct
from django.core.exceptions import ValidationError
from .const import (
BANDTYPE_FLAG_HASNODATA, BANDTYPE_PIXTYPE_MASK, GDAL_TO_POSTGIS,
GDAL_TO_STRUCT, POSTGIS_HEADER_STRUCTURE, POSTGIS_TO_GDAL, STRUCT_SIZE,
)
def pack(structure, data):
"""
Pack data into hex string with little endian format.
"""
return struct.pack('<' + structure, *data)
def unpack(structure, data):
"""
Unpack little endian hexlified binary string into a list.
"""
return struct.unpack('<' + structure, bytes.fromhex(data))
def chunk(data, index):
"""
Split a string into two parts at the input index.
"""
return data[:index], data[index:]
def from_pgraster(data):
"""
Convert a PostGIS HEX String into a dictionary.
"""
if data is None:
return
# Split raster header from data
header, data = chunk(data, 122)
header = unpack(POSTGIS_HEADER_STRUCTURE, header)
# Parse band data
bands = []
pixeltypes = []
while data:
# Get pixel type for this band
pixeltype_with_flags, data = chunk(data, 2)
pixeltype_with_flags = unpack('B', pixeltype_with_flags)[0]
pixeltype = pixeltype_with_flags & BANDTYPE_PIXTYPE_MASK
# Convert datatype from PostGIS to GDAL & get pack type and size
pixeltype = POSTGIS_TO_GDAL[pixeltype]
pack_type = GDAL_TO_STRUCT[pixeltype]
pack_size = 2 * STRUCT_SIZE[pack_type]
# Parse band nodata value. The nodata value is part of the
# PGRaster string even if the nodata flag is True, so it always
# has to be chunked off the data string.
nodata, data = chunk(data, pack_size)
nodata = unpack(pack_type, nodata)[0]
# Chunk and unpack band data (pack size times nr of pixels)
band, data = chunk(data, pack_size * header[10] * header[11])
band_result = {'data': bytes.fromhex(band)}
# Set the nodata value if the nodata flag is set.
if pixeltype_with_flags & BANDTYPE_FLAG_HASNODATA:
band_result['nodata_value'] = nodata
# Append band data to band list
bands.append(band_result)
# Store pixeltype of this band in pixeltypes array
pixeltypes.append(pixeltype)
# Check that all bands have the same pixeltype.
# This is required by GDAL. PostGIS rasters could have different pixeltypes
# for bands of the same raster.
if len(set(pixeltypes)) != 1:
raise ValidationError("Band pixeltypes are not all equal.")
return {
'srid': int(header[9]),
'width': header[10], 'height': header[11],
'datatype': pixeltypes[0],
'origin': (header[5], header[6]),
'scale': (header[3], header[4]),
'skew': (header[7], header[8]),
'bands': bands,
}
def to_pgraster(rast):
"""
Convert a GDALRaster into PostGIS Raster format.
"""
# Prepare the raster header data as a tuple. The first two numbers are
# the endianness and the PostGIS Raster Version, both are fixed by
# PostGIS at the moment.
rasterheader = (
1, 0, len(rast.bands), rast.scale.x, rast.scale.y,
rast.origin.x, rast.origin.y, rast.skew.x, rast.skew.y,
rast.srs.srid, rast.width, rast.height,
)
# Pack raster header.
result = pack(POSTGIS_HEADER_STRUCTURE, rasterheader)
for band in rast.bands:
# The PostGIS raster band header has exactly two elements, a 8BUI byte
# and the nodata value.
#
# The 8BUI stores both the PostGIS pixel data type and a nodata flag.
# It is composed as the datatype with BANDTYPE_FLAG_HASNODATA (1 << 6)
# for existing nodata values:
# 8BUI_VALUE = PG_PIXEL_TYPE (0-11) | BANDTYPE_FLAG_HASNODATA
#
# For example, if the byte value is 71, then the datatype is
# 71 & ~BANDTYPE_FLAG_HASNODATA = 7 (32BSI)
# and the nodata value is True.
structure = 'B' + GDAL_TO_STRUCT[band.datatype()]
# Get band pixel type in PostGIS notation
pixeltype = GDAL_TO_POSTGIS[band.datatype()]
# Set the nodata flag
if band.nodata_value is not None:
pixeltype |= BANDTYPE_FLAG_HASNODATA
# Pack band header
bandheader = pack(structure, (pixeltype, band.nodata_value or 0))
# Add packed header and band data to result
result += bandheader + band.data(as_memoryview=True)
# Convert raster to hex string before passing it to the DB.
return result.hex()

View File

@@ -0,0 +1,71 @@
from django.db.backends.postgresql.schema import DatabaseSchemaEditor
from django.db.models.expressions import Col, Func
class PostGISSchemaEditor(DatabaseSchemaEditor):
geom_index_type = 'GIST'
geom_index_ops_nd = 'GIST_GEOMETRY_OPS_ND'
rast_index_template = 'ST_ConvexHull(%(expressions)s)'
sql_alter_column_to_3d = "ALTER COLUMN %(column)s TYPE %(type)s USING ST_Force3D(%(column)s)::%(type)s"
sql_alter_column_to_2d = "ALTER COLUMN %(column)s TYPE %(type)s USING ST_Force2D(%(column)s)::%(type)s"
def geo_quote_name(self, name):
return self.connection.ops.geo_quote_name(name)
def _field_should_be_indexed(self, model, field):
if getattr(field, 'spatial_index', False):
return True
return super()._field_should_be_indexed(model, field)
def _create_index_sql(self, model, *, fields=None, **kwargs):
if fields is None or len(fields) != 1 or not hasattr(fields[0], 'geodetic'):
return super()._create_index_sql(model, fields=fields, **kwargs)
field = fields[0]
expressions = None
opclasses = None
if field.geom_type == 'RASTER':
# For raster fields, wrap index creation SQL statement with ST_ConvexHull.
# Indexes on raster columns are based on the convex hull of the raster.
expressions = Func(Col(None, field), template=self.rast_index_template)
fields = None
elif field.dim > 2 and not field.geography:
# Use "nd" ops which are fast on multidimensional cases
opclasses = [self.geom_index_ops_nd]
name = kwargs.get('name')
if not name:
name = self._create_index_name(model._meta.db_table, [field.column], '_id')
return super()._create_index_sql(
model,
fields=fields,
name=name,
using=' USING %s' % self.geom_index_type,
opclasses=opclasses,
expressions=expressions,
)
def _alter_column_type_sql(self, table, old_field, new_field, new_type):
"""
Special case when dimension changed.
"""
if not hasattr(old_field, 'dim') or not hasattr(new_field, 'dim'):
return super()._alter_column_type_sql(table, old_field, new_field, new_type)
if old_field.dim == 2 and new_field.dim == 3:
sql_alter = self.sql_alter_column_to_3d
elif old_field.dim == 3 and new_field.dim == 2:
sql_alter = self.sql_alter_column_to_2d
else:
sql_alter = self.sql_alter_column_type
return (
(
sql_alter % {
"column": self.quote_name(new_field.column),
"type": new_type,
},
[],
),
[],
)