test_reproject.py 7.64 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
#!/usr/bin/env python
# -*- coding: utf-8 -*-

"""
test_reproject
--------------

Tests for `py_tools_ds.geo.raster.reproject` module.
"""

import os
from unittest import TestCase

import numpy as np
from gdalnumeric import LoadFile

from py_tools_ds import __path__
from py_tools_ds.geo.coord_calc import corner_coord_to_minmax, get_corner_coordinates
19
from py_tools_ds.geo.raster.reproject import SensorMapGeometryTransformer, SensorMapGeometryTransformer3D
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


tests_path = os.path.abspath(os.path.join(__path__[0], "..", "tests"))


class Test_SensorMapGeometryTransformer(TestCase):
    def setUp(self):
        self.dem_map_geo = LoadFile(os.path.join(tests_path, 'data', 'dem_map_geo.bsq'))
        self.dem_sensor_geo = LoadFile(os.path.join(tests_path, 'data', 'dem_sensor_geo.bsq'))
        self.lons = LoadFile(os.path.join(tests_path, 'data', 'lons_full_vnir.bsq'))
        self.lats = LoadFile(os.path.join(tests_path, 'data', 'lats_full_vnir.bsq'))
        self.dem_area_extent_coarse_subset_utm = [622613.864409047,  # LL_x
                                                  5254111.40255343,  # LL_x
                                                  660473.864409047,  # LL_x
                                                  5269351.40255343]  # UR_y

        self.expected_dem_area_extent_lonlat = [10.685733901515151,  # LL_x
                                                47.44113415492957,  # LL_y
                                                11.073066098484848,  # UR_x
                                                47.54576584507042]  # UR_y

        self.expected_dem_area_extent_utm = [626938.928052,  # LL_x
                                             5256253.56579,  # LL_y
                                             656188.928052,  # UR_x
                                             5267203.56579]  # UR_y

    def test_to_sensor_geometry(self):
47
        SMGT = SensorMapGeometryTransformer(lons=self.lons,
48
49
                                            lats=self.lats,
                                            resamp_alg='nearest')
50
51
        dem_sensors_geo = SMGT.to_sensor_geometry(self.dem_map_geo,
                                                  src_prj=32632, src_extent=self.dem_area_extent_coarse_subset_utm)
52
        self.assertIsInstance(dem_sensors_geo, np.ndarray)
53
        self.assertFalse(np.array_equal(np.unique(dem_sensors_geo), np.array([0])))
54
55
56
        self.assertEquals(dem_sensors_geo.shape, (150, 1000))

    def test_to_map_geometry_lonlat(self):
57
        SMGT = SensorMapGeometryTransformer(lons=self.lons,
58
59
60
61
                                            lats=self.lats,
                                            resamp_alg='nearest')

        # to Lon/Lat
62
        dem_map_geo, dem_gt, dem_prj = SMGT.to_map_geometry(self.dem_sensor_geo, tgt_prj=4326)
63
64
65
66
67
68
69
        self.assertIsInstance(dem_map_geo, np.ndarray)
        self.assertEquals(dem_map_geo.shape, (286, 1058))
        xmin, xmax, ymin, ymax = corner_coord_to_minmax(get_corner_coordinates(gt=dem_gt,
                                                                               cols=dem_map_geo.shape[1],
                                                                               rows=dem_map_geo.shape[0]))
        self.assertTrue(False not in np.isclose(np.array([xmin, ymin, xmax, ymax]),
                                                np.array(self.expected_dem_area_extent_lonlat)))
70
        self.assertFalse(np.array_equal(np.unique(dem_map_geo), np.array([0])))
71

72
73
74
        with self.assertRaises(ValueError):
            SMGT.to_map_geometry(self.dem_sensor_geo[:10, :10], tgt_prj=4326)  # must have the shape of lons/lats

75
    def test_to_map_geometry_utm(self):
76
        SMGT = SensorMapGeometryTransformer(lons=self.lons,
77
78
79
80
                                            lats=self.lats,
                                            resamp_alg='nearest')

        # to UTM32
81
        dem_map_geo, dem_gt, dem_prj = SMGT.to_map_geometry(self.dem_sensor_geo, tgt_prj=32632, tgt_res=(30, 30))
82
83
84
85
86
87
88
        self.assertIsInstance(dem_map_geo, np.ndarray)
        self.assertEquals(dem_map_geo.shape, (365, 975))
        xmin, xmax, ymin, ymax = corner_coord_to_minmax(get_corner_coordinates(gt=dem_gt,
                                                                               cols=dem_map_geo.shape[1],
                                                                               rows=dem_map_geo.shape[0]))
        self.assertTrue(False not in np.isclose(np.array([xmin, ymin, xmax, ymax]),
                                                np.array(self.expected_dem_area_extent_utm)))
89
        self.assertFalse(np.array_equal(np.unique(dem_map_geo), np.array([0])))
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


class Test_SensorMapGeometryTransformer3D(TestCase):
    def setUp(self):
        dem_map_geo = LoadFile(os.path.join(tests_path, 'data', 'dem_map_geo.bsq'))
        dem_sensor_geo = LoadFile(os.path.join(tests_path, 'data', 'dem_sensor_geo.bsq'))
        lons = LoadFile(os.path.join(tests_path, 'data', 'lons_full_vnir.bsq'))
        lats = LoadFile(os.path.join(tests_path, 'data', 'lats_full_vnir.bsq'))

        self.data_map_geo_3D = np.dstack([dem_map_geo, dem_map_geo])
        self.data_sensor_geo_3D = np.dstack([dem_sensor_geo, dem_sensor_geo])
        self.lons_3D = np.dstack([lons, lons])  # TODO use different lons per band here
        self.lats_3D = np.dstack([lats, lats])  # TODO use different lats per band here

        self.dem_area_extent_coarse_subset_utm = [622613.864409047,  # LL_x
                                                  5254111.40255343,  # LL_x
                                                  660473.864409047,  # LL_x
                                                  5269351.40255343]  # UR_y

        self.expected_dem_area_extent_lonlat = [10.685733901515151,  # LL_x
                                                47.44113415492957,  # LL_y
                                                11.073066098484848,  # UR_x
                                                47.54576584507042]  # UR_y

        self.expected_dem_area_extent_utm = [626938.928052,  # LL_x
                                             5256253.56579,  # LL_y
                                             656188.928052,  # UR_x
                                             5267203.56579]  # UR_y

    def test_to_map_geometry_lonlat_3D_geolayer(self):
        SMGT = SensorMapGeometryTransformer3D(lons=self.lons_3D,
                                              lats=self.lats_3D,
                                              resamp_alg='nearest')

        # to Lon/Lat
        data_mapgeo_3D, dem_gt, dem_prj = SMGT.to_map_geometry(self.data_sensor_geo_3D, tgt_prj=4326)
        self.assertIsInstance(data_mapgeo_3D, np.ndarray)
        self.assertEquals(data_mapgeo_3D.shape, (286, 1058, 2))
        xmin, xmax, ymin, ymax = corner_coord_to_minmax(get_corner_coordinates(gt=dem_gt,
                                                                               cols=data_mapgeo_3D.shape[1],
                                                                               rows=data_mapgeo_3D.shape[0]))
        self.assertTrue(False not in np.isclose(np.array([xmin, ymin, xmax, ymax]),
                                                np.array(self.expected_dem_area_extent_lonlat)))

    def test_to_sensor_geometry(self):
        SMGT = SensorMapGeometryTransformer3D(lons=self.lons_3D,
                                              lats=self.lats_3D,
                                              resamp_alg='nearest')
        dem_sensors_geo = SMGT.to_sensor_geometry(self.data_map_geo_3D,
                                                  src_prj=32632, src_extent=self.dem_area_extent_coarse_subset_utm)
        self.assertIsInstance(dem_sensors_geo, np.ndarray)
        self.assertEquals(dem_sensors_geo.shape, (150, 1000,  2))