1"""
2Reproject a Raster using ST_Transform
3=====================================
4
5The `ST_Transform()` function (and a few others like `ST_SnapToGrid()`) can be used on
6both `Geometry` and `Raster` types. In `GeoAlchemy2`, this function is only defined for
7`Geometry` as it can not be defined for several types at the same time. Thus using this
8function on `Raster` requires minor tweaking.
9
10This example uses both SQLAlchemy core and ORM queries.
11"""
12from sqlalchemy import Column
13from sqlalchemy import func
14from sqlalchemy import Integer
15from sqlalchemy import MetaData
16from sqlalchemy import select
17from sqlalchemy import Table
18from sqlalchemy.orm import Query
19from sqlalchemy.ext.declarative import declarative_base
20
21from geoalchemy2 import Geometry
22from geoalchemy2 import Raster
23
24
25metadata = MetaData()
26Base = declarative_base(metadata=metadata)
27
28table = Table(
29    "raster_table",
30    metadata,
31    Column("id", Integer, primary_key=True),
32    Column("geom", Geometry("POLYGON", 4326)),
33    Column("rast", Raster(srid=4326)),
34)
35
36
37class RasterTable(Base):
38    __tablename__ = 'raster_table_orm'
39    id = Column(Integer, primary_key=True)
40    geom = Column(Geometry("POLYGON", 4326))
41    rast = Column(Raster(srid=4326))
42
43    def __init__(self, rast):
44        self.rast = rast
45
46
47def test_transform_core():
48    # Define the transform query for both the geometry and the raster in a naive way
49    wrong_query = select([
50        func.ST_Transform(table.c.geom, 2154),
51        func.ST_Transform(table.c.rast, 2154)
52    ])
53
54    # Check the query
55    assert str(wrong_query) == (
56        "SELECT "
57        "ST_AsEWKB("
58        "ST_Transform(raster_table.geom, :ST_Transform_2)) AS \"ST_Transform_1\", "
59        "ST_AsEWKB("  # <= Note that the raster is processed as a Geometry here
60        "ST_Transform(raster_table.rast, :ST_Transform_4)) AS \"ST_Transform_3\" \n"
61        "FROM raster_table"
62    )
63
64    # Define the transform query for both the geometry and the raster in the correct way
65    correct_query = select([
66        func.ST_Transform(table.c.geom, 2154),
67        func.ST_Transform(table.c.rast, 2154, type_=Raster)
68    ])
69
70    # Check the query
71    assert str(correct_query) == (
72        "SELECT "
73        "ST_AsEWKB("
74        "ST_Transform(raster_table.geom, :ST_Transform_2)) AS \"ST_Transform_1\", "
75        "raster("  # <= This time the raster is correctly processed as a Raster
76        "ST_Transform(raster_table.rast, :ST_Transform_4)) AS \"ST_Transform_3\" \n"
77        "FROM raster_table"
78    )
79
80
81def test_transform_ORM():
82    # Define the transform query for both the geometry and the raster in a naive way
83    wrong_query = Query([
84        RasterTable.geom.ST_Transform(2154),
85        RasterTable.rast.ST_Transform(2154)
86    ])
87
88    # Check the query
89    assert str(wrong_query) == (
90        "SELECT "
91        "ST_AsEWKB("
92        "ST_Transform(raster_table_orm.geom, :ST_Transform_2)) AS \"ST_Transform_1\", "
93        "ST_AsEWKB("  # <= Note that the raster is processed as a Geometry here
94        "ST_Transform(raster_table_orm.rast, :ST_Transform_4)) AS \"ST_Transform_3\" \n"
95        "FROM raster_table_orm"
96    )
97
98    # Define the transform query for both the geometry and the raster in the correct way
99    correct_query = Query([
100        RasterTable.geom.ST_Transform(2154),
101        RasterTable.rast.ST_Transform(2154, type_=Raster)
102    ])
103
104    # Check the query
105    assert str(correct_query) == (
106        "SELECT "
107        "ST_AsEWKB("
108        "ST_Transform(raster_table_orm.geom, :ST_Transform_2)) AS \"ST_Transform_1\", "
109        "raster("  # <= This time the raster is correctly processed as a Raster
110        "ST_Transform(raster_table_orm.rast, :ST_Transform_4)) AS \"ST_Transform_3\" \n"
111        "FROM raster_table_orm"
112    )
113