tes.ipynb
{
"cells": [
{
"cell_type": "code",
"execution_count": 31,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import rasterio\n",
"from typing import Tuple, Union, List, Optional, Any\n",
"from affine import Affine\n",
"from rasterio.crs import CRS\n",
"from rasterio.plot import reshape_as_image"
]
},
{
"cell_type": "code",
"execution_count": 49,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": "array([[50, 51, 52, 53, 54, 55, 56, 57, 58, 59],\n [60, 61, 62, 63, 64, 65, 66, 67, 68, 69],\n [70, 71, 72, 73, 74, 75, 76, 77, 78, 79]])"
},
"execution_count": 49,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"np.arange(10*10).reshape(10,10)[-5:8]"
]
},
{
"cell_type": "code",
"execution_count": 51,
"metadata": {},
"outputs": [],
"source": [
"def parse_slicer(key: Union[int, slice, None], length:int) -> int:\n",
" if key is None:\n",
" start = 0\n",
" elif isinstance(key, int):\n",
" start = key if key > 0 else length + key\n",
" elif isinstance(key, slice):\n",
" if slice.start is None:\n",
" start = 0\n",
" elif slice.start < 0:\n",
" start = length + slice.start\n",
" elif slice.start > 0:\n",
" start = slice.start\n",
" else:\n",
" raise ValueError\n",
" else:\n",
" raise ValueError\n",
" return start"
]
},
{
"cell_type": "code",
"execution_count": 52,
"metadata": {},
"outputs": [],
"source": [
"class A(np.ndarray):\n",
" def __new__(cls, array: np.ndarray):\n",
" return array.view(cls)\n",
" def __getitem__(self, keys):\n",
" if len(keys) == 1:\n",
" row_col_min: List[int] = [parse_slicer(keys, self.__array__().shape[0]),0]\n",
" elif len(keys) == 2:\n",
" row_col_min = [parse_slicer(key, self.__array__().shape[i]) for i,key in enumerate(keys)]\n",
" elif len(keys) == 3:\n",
" row_col_min = [parse_slicer(key, self.__array__().shape[i]) for i,key in enumerate(keys[:2])]\n",
"\n",
" print(row_col_min)\n",
" \n",
" return self.__array__()[keys]"
]
},
{
"cell_type": "code",
"execution_count": 43,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[9, 1]\n"
]
},
{
"data": {
"text/plain": "array([], dtype=int64)"
},
"execution_count": 43,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"A(np.arange(10*10).reshape(10,10))[-1:4,1]"
]
},
{
"cell_type": "code",
"execution_count": 44,
"metadata": {},
"outputs": [
{
"ename": "SyntaxError",
"evalue": "invalid syntax (<ipython-input-44-f81480c15b4a>, line 1)",
"output_type": "error",
"traceback": [
"\u001b[0;36m File \u001b[0;32m\"<ipython-input-44-f81480c15b4a>\"\u001b[0;36m, line \u001b[0;32m1\u001b[0m\n\u001b[0;31m -1:\u001b[0m\n\u001b[0m ^\u001b[0m\n\u001b[0;31mSyntaxError\u001b[0m\u001b[0;31m:\u001b[0m invalid syntax\n"
]
}
],
"source": [
"-1:"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {},
"outputs": [],
"source": [
"class Raster(np.ndarray):\n",
" def __init__(\n",
" self,\n",
" array: np.ndarray,\n",
" resolution: Union[\n",
" None, Tuple[float, float], List[float], Tuple[float, ...], float\n",
" ] = None,\n",
" x_min: Optional[float] = None,\n",
" y_max: Optional[float] = None,\n",
" x_max: Optional[float] = None,\n",
" y_min: Optional[float] = None,\n",
" epsg: int = 4326,\n",
" no_data: Union[float, int] = -32767.0,\n",
" transform: Optional[Affine] = None,\n",
" source: Optional[str] = None,\n",
" ):\n",
" if transform is None:\n",
" if resolution is None and x_min is None and y_min is None:\n",
" raise ValueError(\n",
" \"Please define resolution and at least x minimum and y minimum\"\n",
" )\n",
"\n",
" if resolution is not None and x_min is None and y_max is None:\n",
" raise ValueError(\"Please define x_min and y_max\")\n",
"\n",
" if isinstance(resolution, float):\n",
" self.resolution: Tuple[float, float] = (\n",
" resolution,\n",
" resolution,\n",
" )\n",
" elif isinstance(resolution, Iterable):\n",
" self.resolution = (resolution[0], resolution[1])\n",
"\n",
" if (\n",
" resolution is None\n",
" and x_min is not None\n",
" and y_min is not None\n",
" and x_max is not None\n",
" and y_max is not None\n",
" ):\n",
" self.resolution = (\n",
" (x_max - x_min) / array.shape[1],\n",
" (y_max - y_min) / array.shape[0],\n",
" )\n",
"\n",
" self.transform: Affine = Affine.translation(x_min, y_max) * Affine.scale(\n",
" self.resolution[0], -self.resolution[1]\n",
" )\n",
" elif isinstance(transform, Affine):\n",
" self.transform = transform\n",
" self.resolution = (transform[0], abs(transform[4]))\n",
" else:\n",
" raise ValueError(\n",
" \"Please define affine parameter or resolution and xmin ymax\"\n",
" )\n",
"\n",
" self.epsg = epsg\n",
"\n",
" self.crs = CRS.from_epsg(epsg)\n",
" self.no_data = no_data\n",
" self.source = source\n",
" self.__check_validity()\n",
"\n",
" def __new__(cls, array: np.ndarray, *args, **kwargs) -> \"Raster\":\n",
" return array.view(cls)\n",
"\n",
" def __getitem__(self, key: Union[int, Tuple[Any, ...], slice]) -> np.ndarray:\n",
" if len(keys) == 1:\n",
" row_col_min: List[int] = [parse_slicer(keys),0]\n",
" elif len(keys) == 2:\n",
" row_col_min = [parse_slicer(key) for key in keys]\n",
" elif len(keys) == 3:\n",
" row_col_min = [parse_slicer(key) for key in keys[:2]]\n",
" return self.array.__getitem__(key)\n",
"\n",
" @classmethod\n",
" def from_rasterfile(cls, raster_file: str) -> \"Raster\":\n",
" \"\"\"Get raster from supported gdal raster file\n",
"\n",
" Parameters\n",
" -------\n",
" raster_file : str\n",
" location of raser file\n",
"\n",
" Returns\n",
" -------\n",
" Raster\n",
" \"\"\"\n",
" with rasterio.open(raster_file) as file:\n",
" _raster = reshape_as_image(file.read())\n",
"\n",
" return cls(\n",
" _raster,\n",
" transform=file.transform,\n",
" epsg=file.crs.to_epsg(),\n",
" no_data=file.nodatavals[0],\n",
" source=raster_file,\n",
" )\n",
"\n",
" @property\n",
" def array(self) -> np.ndarray:\n",
" \"\"\"the numpy array of raster\"\"\"\n",
" return self.__array__()\n",
"\n",
" @property\n",
" def __transform(self) -> Tuple[float, ...]:\n",
" return tuple(self.transform)\n",
"\n",
" @property\n",
" def x_min(self) -> float:\n",
" \"\"\"minimum x-axis coordinate\"\"\"\n",
" return self.__transform[2]\n",
"\n",
" @property\n",
" def y_max(self) -> float:\n",
" \"\"\"maximum y-axis coordinate\"\"\"\n",
" return self.__transform[5]\n",
"\n",
" @property\n",
" def x_max(self) -> float:\n",
" \"\"\"maximum x-axis coordinate\"\"\"\n",
" return self.__transform[2] + (self.resolution[0] * self.cols)\n",
"\n",
" @property\n",
" def y_min(self) -> float:\n",
" \"\"\"minimum y-axis coordinate\"\"\"\n",
" return self.__transform[5] - (self.resolution[1] * self.rows)\n",
"\n",
" @property\n",
" def top(self) -> float:\n",
" \"\"\"top y-axis coordinate\"\"\"\n",
" return self.y_max\n",
"\n",
" @property\n",
" def left(self) -> float:\n",
" \"\"\"left x-axis coordinate\"\"\"\n",
" return self.x_min\n",
"\n",
" @property\n",
" def right(self) -> float:\n",
" \"\"\"right x-axis coordinate\"\"\"\n",
" return self.x_max\n",
"\n",
" @property\n",
" def bottom(self) -> float:\n",
" \"\"\"bottom y-axis coordinate\"\"\"\n",
" return self.y_min\n",
"\n",
" @property\n",
" def rows(self) -> int:\n",
" \"\"\"number of row, height\"\"\"\n",
" return int(self.array.shape[0])\n",
"\n",
" @property\n",
" def cols(self) -> int:\n",
" \"\"\"number of column, width\"\"\"\n",
" return int(self.array.shape[1])\n",
"\n",
" @property\n",
" def layers(self) -> int:\n",
" \"\"\"number of layer / channel\"\"\"\n",
" _layers: int = 1\n",
" if len(self.array.shape) > 2:\n",
" _layers = self.array.shape[2]\n",
" return _layers\n",
"\n",
" @property\n",
" def x_extent(self) -> float:\n",
" \"\"\"width of raster in the map unit (degree decimal or meters)\"\"\"\n",
" return self.x_max - self.x_min\n",
"\n",
" @property\n",
" def y_extent(self) -> float:\n",
" \"\"\"height of raster in the map unit (degree decimal or meters)\"\"\"\n",
" return self.y_max - self.y_min\n",
"\n",
" @property\n",
" def is_projected(self) -> bool:\n",
" \"\"\"check crs is projected or not\"\"\"\n",
" return self.crs.is_projected\n",
"\n",
" @property\n",
" def is_geographic(self) -> bool:\n",
" \"\"\"check crs is geographic or not\"\"\"\n",
" return self.crs.is_geographic\n",
"\n",
" def __check_validity(self) -> None:\n",
" \"\"\"Check geometry validity\n",
"\n",
" Raises\n",
" ------\n",
" ValueError\n",
" x min, y min is greater than x max, y max\n",
" ValueError\n",
" x min is greater than x max\n",
" ValueError\n",
" y min is greater than y max\n",
" \"\"\"\n",
" if self.x_extent < 0 and self.y_extent < 0:\n",
" raise ValueError(\n",
" \"x min should be less than x max and y min should be less than y max\"\n",
" )\n",
" elif self.x_extent < 0 and self.y_extent > 0:\n",
" raise ValueError(\"x min should be less than x max\")\n",
" elif self.x_extent > 0 and self.y_extent < 0:\n",
" raise ValueError(\"y min should be less than y max\")"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": "(8668, 2721, 1)"
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Raster.from_rasterfile(\"/mnt/d/temp/ncku.tif\").shape"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3.8.5 64-bit ('.venv': venv)",
"metadata": {
"interpreter": {
"hash": "b2a74825ea9f7432a70215a3bbf56341fecb286da4b95abaafab5c59fae8d59d"
}
},
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.5"
},
"orig_nbformat": 2
},
"nbformat": 4,
"nbformat_minor": 2
}