You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

pgraster.py 4.3KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141
  1. import struct
  2. from django.forms import ValidationError
  3. from .const import (
  4. GDAL_TO_POSTGIS, GDAL_TO_STRUCT, POSTGIS_HEADER_STRUCTURE, POSTGIS_TO_GDAL,
  5. STRUCT_SIZE,
  6. )
  7. def pack(structure, data):
  8. """
  9. Pack data into hex string with little endian format.
  10. """
  11. return struct.pack('<' + structure, *data)
  12. def unpack(structure, data):
  13. """
  14. Unpack little endian hexlified binary string into a list.
  15. """
  16. return struct.unpack('<' + structure, bytes.fromhex(data))
  17. def chunk(data, index):
  18. """
  19. Split a string into two parts at the input index.
  20. """
  21. return data[:index], data[index:]
  22. def from_pgraster(data):
  23. """
  24. Convert a PostGIS HEX String into a dictionary.
  25. """
  26. if data is None:
  27. return
  28. # Split raster header from data
  29. header, data = chunk(data, 122)
  30. header = unpack(POSTGIS_HEADER_STRUCTURE, header)
  31. # Parse band data
  32. bands = []
  33. pixeltypes = []
  34. while data:
  35. # Get pixel type for this band
  36. pixeltype, data = chunk(data, 2)
  37. pixeltype = unpack('B', pixeltype)[0]
  38. # Subtract nodata byte from band nodata value if it exists
  39. has_nodata = pixeltype >= 64
  40. if has_nodata:
  41. pixeltype -= 64
  42. # Convert datatype from PostGIS to GDAL & get pack type and size
  43. pixeltype = POSTGIS_TO_GDAL[pixeltype]
  44. pack_type = GDAL_TO_STRUCT[pixeltype]
  45. pack_size = 2 * STRUCT_SIZE[pack_type]
  46. # Parse band nodata value. The nodata value is part of the
  47. # PGRaster string even if the nodata flag is True, so it always
  48. # has to be chunked off the data string.
  49. nodata, data = chunk(data, pack_size)
  50. nodata = unpack(pack_type, nodata)[0]
  51. # Chunk and unpack band data (pack size times nr of pixels)
  52. band, data = chunk(data, pack_size * header[10] * header[11])
  53. band_result = {'data': bytes.fromhex(band)}
  54. # If the nodata flag is True, set the nodata value.
  55. if has_nodata:
  56. band_result['nodata_value'] = nodata
  57. # Append band data to band list
  58. bands.append(band_result)
  59. # Store pixeltype of this band in pixeltypes array
  60. pixeltypes.append(pixeltype)
  61. # Check that all bands have the same pixeltype.
  62. # This is required by GDAL. PostGIS rasters could have different pixeltypes
  63. # for bands of the same raster.
  64. if len(set(pixeltypes)) != 1:
  65. raise ValidationError("Band pixeltypes are not all equal.")
  66. return {
  67. 'srid': int(header[9]),
  68. 'width': header[10], 'height': header[11],
  69. 'datatype': pixeltypes[0],
  70. 'origin': (header[5], header[6]),
  71. 'scale': (header[3], header[4]),
  72. 'skew': (header[7], header[8]),
  73. 'bands': bands,
  74. }
  75. def to_pgraster(rast):
  76. """
  77. Convert a GDALRaster into PostGIS Raster format.
  78. """
  79. # Prepare the raster header data as a tuple. The first two numbers are
  80. # the endianness and the PostGIS Raster Version, both are fixed by
  81. # PostGIS at the moment.
  82. rasterheader = (
  83. 1, 0, len(rast.bands), rast.scale.x, rast.scale.y,
  84. rast.origin.x, rast.origin.y, rast.skew.x, rast.skew.y,
  85. rast.srs.srid, rast.width, rast.height,
  86. )
  87. # Pack raster header.
  88. result = pack(POSTGIS_HEADER_STRUCTURE, rasterheader)
  89. for band in rast.bands:
  90. # The PostGIS raster band header has exactly two elements, a 8BUI byte
  91. # and the nodata value.
  92. #
  93. # The 8BUI stores both the PostGIS pixel data type and a nodata flag.
  94. # It is composed as the datatype integer plus 64 as a flag for existing
  95. # nodata values:
  96. # 8BUI_VALUE = PG_PIXEL_TYPE (0-11) + FLAG (0 or 64)
  97. #
  98. # For example, if the byte value is 71, then the datatype is
  99. # 71-64 = 7 (32BSI) and the nodata value is True.
  100. structure = 'B' + GDAL_TO_STRUCT[band.datatype()]
  101. # Get band pixel type in PostGIS notation
  102. pixeltype = GDAL_TO_POSTGIS[band.datatype()]
  103. # Set the nodata flag
  104. if band.nodata_value is not None:
  105. pixeltype += 64
  106. # Pack band header
  107. bandheader = pack(structure, (pixeltype, band.nodata_value or 0))
  108. # Add packed header and band data to result
  109. result += bandheader + band.data(as_memoryview=True)
  110. # Convert raster to hex string before passing it to the DB.
  111. return result.hex()