band.py 8.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273
  1. from ctypes import byref, c_double, c_int, c_void_p
  2. from django.contrib.gis.gdal.error import GDALException
  3. from django.contrib.gis.gdal.prototypes import raster as capi
  4. from django.contrib.gis.gdal.raster.base import GDALRasterBase
  5. from django.contrib.gis.shortcuts import numpy
  6. from django.utils.encoding import force_str
  7. from .const import (
  8. GDAL_COLOR_TYPES,
  9. GDAL_INTEGER_TYPES,
  10. GDAL_PIXEL_TYPES,
  11. GDAL_TO_CTYPES,
  12. )
  13. class GDALBand(GDALRasterBase):
  14. """
  15. Wrap a GDAL raster band, needs to be obtained from a GDALRaster object.
  16. """
  17. def __init__(self, source, index):
  18. self.source = source
  19. self._ptr = capi.get_ds_raster_band(source._ptr, index)
  20. def _flush(self):
  21. """
  22. Call the flush method on the Band's parent raster and force a refresh
  23. of the statistics attribute when requested the next time.
  24. """
  25. self.source._flush()
  26. self._stats_refresh = True
  27. @property
  28. def description(self):
  29. """
  30. Return the description string of the band.
  31. """
  32. return force_str(capi.get_band_description(self._ptr))
  33. @property
  34. def width(self):
  35. """
  36. Width (X axis) in pixels of the band.
  37. """
  38. return capi.get_band_xsize(self._ptr)
  39. @property
  40. def height(self):
  41. """
  42. Height (Y axis) in pixels of the band.
  43. """
  44. return capi.get_band_ysize(self._ptr)
  45. @property
  46. def pixel_count(self):
  47. """
  48. Return the total number of pixels in this band.
  49. """
  50. return self.width * self.height
  51. _stats_refresh = False
  52. def statistics(self, refresh=False, approximate=False):
  53. """
  54. Compute statistics on the pixel values of this band.
  55. The return value is a tuple with the following structure:
  56. (minimum, maximum, mean, standard deviation).
  57. If approximate=True, the statistics may be computed based on overviews
  58. or a subset of image tiles.
  59. If refresh=True, the statistics will be computed from the data directly,
  60. and the cache will be updated where applicable.
  61. For empty bands (where all pixel values are nodata), all statistics
  62. values are returned as None.
  63. For raster formats using Persistent Auxiliary Metadata (PAM) services,
  64. the statistics might be cached in an auxiliary file.
  65. """
  66. # Prepare array with arguments for capi function
  67. smin, smax, smean, sstd = c_double(), c_double(), c_double(), c_double()
  68. stats_args = [
  69. self._ptr,
  70. c_int(approximate),
  71. byref(smin),
  72. byref(smax),
  73. byref(smean),
  74. byref(sstd),
  75. c_void_p(),
  76. c_void_p(),
  77. ]
  78. if refresh or self._stats_refresh:
  79. func = capi.compute_band_statistics
  80. else:
  81. # Add additional argument to force computation if there is no
  82. # existing PAM file to take the values from.
  83. force = True
  84. stats_args.insert(2, c_int(force))
  85. func = capi.get_band_statistics
  86. # Computation of statistics fails for empty bands.
  87. try:
  88. func(*stats_args)
  89. result = smin.value, smax.value, smean.value, sstd.value
  90. except GDALException:
  91. result = (None, None, None, None)
  92. self._stats_refresh = False
  93. return result
  94. @property
  95. def min(self):
  96. """
  97. Return the minimum pixel value for this band.
  98. """
  99. return self.statistics()[0]
  100. @property
  101. def max(self):
  102. """
  103. Return the maximum pixel value for this band.
  104. """
  105. return self.statistics()[1]
  106. @property
  107. def mean(self):
  108. """
  109. Return the mean of all pixel values of this band.
  110. """
  111. return self.statistics()[2]
  112. @property
  113. def std(self):
  114. """
  115. Return the standard deviation of all pixel values of this band.
  116. """
  117. return self.statistics()[3]
  118. @property
  119. def nodata_value(self):
  120. """
  121. Return the nodata value for this band, or None if it isn't set.
  122. """
  123. # Get value and nodata exists flag
  124. nodata_exists = c_int()
  125. value = capi.get_band_nodata_value(self._ptr, nodata_exists)
  126. if not nodata_exists:
  127. value = None
  128. # If the pixeltype is an integer, convert to int
  129. elif self.datatype() in GDAL_INTEGER_TYPES:
  130. value = int(value)
  131. return value
  132. @nodata_value.setter
  133. def nodata_value(self, value):
  134. """
  135. Set the nodata value for this band.
  136. """
  137. if value is None:
  138. capi.delete_band_nodata_value(self._ptr)
  139. elif not isinstance(value, (int, float)):
  140. raise ValueError("Nodata value must be numeric or None.")
  141. else:
  142. capi.set_band_nodata_value(self._ptr, value)
  143. self._flush()
  144. def datatype(self, as_string=False):
  145. """
  146. Return the GDAL Pixel Datatype for this band.
  147. """
  148. dtype = capi.get_band_datatype(self._ptr)
  149. if as_string:
  150. dtype = GDAL_PIXEL_TYPES[dtype]
  151. return dtype
  152. def color_interp(self, as_string=False):
  153. """Return the GDAL color interpretation for this band."""
  154. color = capi.get_band_color_interp(self._ptr)
  155. if as_string:
  156. color = GDAL_COLOR_TYPES[color]
  157. return color
  158. def data(self, data=None, offset=None, size=None, shape=None, as_memoryview=False):
  159. """
  160. Read or writes pixel values for this band. Blocks of data can
  161. be accessed by specifying the width, height and offset of the
  162. desired block. The same specification can be used to update
  163. parts of a raster by providing an array of values.
  164. Allowed input data types are bytes, memoryview, list, tuple, and array.
  165. """
  166. offset = offset or (0, 0)
  167. size = size or (self.width - offset[0], self.height - offset[1])
  168. shape = shape or size
  169. if any(x <= 0 for x in size):
  170. raise ValueError("Offset too big for this raster.")
  171. if size[0] > self.width or size[1] > self.height:
  172. raise ValueError("Size is larger than raster.")
  173. # Create ctypes type array generator
  174. ctypes_array = GDAL_TO_CTYPES[self.datatype()] * (shape[0] * shape[1])
  175. if data is None:
  176. # Set read mode
  177. access_flag = 0
  178. # Prepare empty ctypes array
  179. data_array = ctypes_array()
  180. else:
  181. # Set write mode
  182. access_flag = 1
  183. # Instantiate ctypes array holding the input data
  184. if isinstance(data, (bytes, memoryview)) or (
  185. numpy and isinstance(data, numpy.ndarray)
  186. ):
  187. data_array = ctypes_array.from_buffer_copy(data)
  188. else:
  189. data_array = ctypes_array(*data)
  190. # Access band
  191. capi.band_io(
  192. self._ptr,
  193. access_flag,
  194. offset[0],
  195. offset[1],
  196. size[0],
  197. size[1],
  198. byref(data_array),
  199. shape[0],
  200. shape[1],
  201. self.datatype(),
  202. 0,
  203. 0,
  204. )
  205. # Return data as numpy array if possible, otherwise as list
  206. if data is None:
  207. if as_memoryview:
  208. return memoryview(data_array)
  209. elif numpy:
  210. # reshape() needs a reshape parameter with the height first.
  211. return numpy.frombuffer(
  212. data_array, dtype=numpy.dtype(data_array)
  213. ).reshape(tuple(reversed(size)))
  214. else:
  215. return list(data_array)
  216. else:
  217. self._flush()
  218. class BandList(list):
  219. def __init__(self, source):
  220. self.source = source
  221. super().__init__()
  222. def __iter__(self):
  223. for idx in range(1, len(self) + 1):
  224. yield GDALBand(self.source, idx)
  225. def __len__(self):
  226. return capi.get_ds_raster_count(self.source._ptr)
  227. def __getitem__(self, index):
  228. try:
  229. return GDALBand(self.source, index + 1)
  230. except GDALException:
  231. raise GDALException("Unable to get band index %d" % index)