]> granicus.if.org Git - postgis/commitdiff
Removed requirements of SRID for calling GDAL Warp API
authorBborie Park <bkpark at ucdavis.edu>
Sat, 1 Dec 2012 01:19:45 +0000 (01:19 +0000)
committerBborie Park <bkpark at ucdavis.edu>
Sat, 1 Dec 2012 01:19:45 +0000 (01:19 +0000)
git-svn-id: http://svn.osgeo.org/postgis/trunk@10773 b70326c6-7e19-0410-871a-916f4a2858ee

raster/rt_core/rt_api.c
raster/rt_core/rt_api.h
raster/rt_pg/rt_pg.c
raster/test/regress/rt_gdalwarp.sql
raster/test/regress/rt_gdalwarp_expected

index 763089424b99a9c72de90134c32fe6929c0d8e5e..07f1be1cb02dcfb3b6a4a0b2850e6cc92ebab050 100644 (file)
@@ -9240,8 +9240,8 @@ rt_raster_from_gdal_dataset(GDALDatasetH ds) {
  * @return the warped raster or NULL
  */
 rt_raster rt_raster_gdal_warp(
-       rt_raster raster, const char *src_srs,
-       const char *dst_srs,
+       rt_raster raster,
+       const char *src_srs, const char *dst_srs,
        double *scale_x, double *scale_y,
        int *width, int *height,
        double *ul_xw, double *ul_yw,
@@ -9288,7 +9288,6 @@ rt_raster rt_raster_gdal_warp(
        RASTER_DEBUG(3, "starting");
 
        assert(NULL != raster);
-       assert(NULL != src_srs);
 
        /*
                max_err must be gte zero
@@ -9298,20 +9297,34 @@ rt_raster rt_raster_gdal_warp(
        if (max_err < 0.) max_err = 0.125;
        RASTER_DEBUGF(4, "max_err = %f", max_err);
 
-       _src_srs = rt_util_gdal_convert_sr(src_srs, 0);
-       /* dst_srs not provided, set to src_srs */
-       if (NULL == dst_srs)
-               _dst_srs = rt_util_gdal_convert_sr(src_srs, 0);
-       else
-               _dst_srs = rt_util_gdal_convert_sr(dst_srs, 0);
+       /* handle srs */
+       if (src_srs != NULL) {
+               /* reprojection taking place */
+               if (dst_srs != NULL && strcmp(src_srs, dst_srs) != 0) {
+                       RASTER_DEBUG(4, "Warp operation does include a reprojection");
+                       _src_srs = rt_util_gdal_convert_sr(src_srs, 0);
+                       _dst_srs = rt_util_gdal_convert_sr(dst_srs, 0);
+               }
+               /* no reprojection, a stub just for clarity */
+               else {
+                       RASTER_DEBUG(4, "Warp operation does NOT include reprojection");
+               }
+       }
+       else if (dst_srs != NULL) {
+               /* dst_srs provided but not src_srs */
+               rterror("rt_raster_gdal_warp: SRS required for input raster if SRS provided for warped raster");
+               return NULL;
+       }
 
        /* load raster into a GDAL MEM dataset */
        src_ds = rt_raster_to_gdal_mem(raster, _src_srs, NULL, NULL, 0, &src_drv);
        if (NULL == src_ds) {
                rterror("rt_raster_gdal_warp: Unable to convert raster to GDAL MEM format");
 
-               CPLFree(_src_srs);
-               CPLFree(_dst_srs);
+               if (_src_srs != NULL) {
+                       CPLFree(_src_srs);
+                       CPLFree(_dst_srs);
+               }
 
                return NULL;
        }
@@ -9324,18 +9337,26 @@ rt_raster rt_raster_gdal_warp(
 
                GDALClose(src_ds);
 
-               CPLFree(_src_srs);
-               CPLFree(_dst_srs);
+               if (_src_srs != NULL) {
+                       CPLFree(_src_srs);
+                       CPLFree(_dst_srs);
+               }
 
                return NULL;
        }
        for (i = 0; i < transform_opts_len; i++) {
                switch (i) {
                        case 1:
-                               transform_opts[i] = (char *) rtalloc(sizeof(char) * (strlen("DST_SRS=") + strlen(_dst_srs) + 1));
+                               if (_dst_srs != NULL)
+                                       transform_opts[i] = (char *) rtalloc(sizeof(char) * (strlen("DST_SRS=") + strlen(_dst_srs) + 1));
+                               else
+                                       transform_opts[i] = (char *) rtalloc(sizeof(char) * (strlen("DST_SRS=") + 1));
                                break;
                        case 0:
-                               transform_opts[i] = (char *) rtalloc(sizeof(char) * (strlen("SRC_SRS=") + strlen(_src_srs) + 1));
+                               if (_src_srs != NULL)
+                                       transform_opts[i] = (char *) rtalloc(sizeof(char) * (strlen("SRC_SRS=") + strlen(_src_srs) + 1));
+                               else
+                                       transform_opts[i] = (char *) rtalloc(sizeof(char) * (strlen("SRC_SRS=") + 1));
                                break;
                }
                if (NULL == transform_opts[i]) {
@@ -9345,35 +9366,44 @@ rt_raster rt_raster_gdal_warp(
                        rtdealloc(transform_opts);
                        GDALClose(src_ds);
 
-                       CPLFree(_src_srs);
-                       CPLFree(_dst_srs);
+                       if (_src_srs != NULL) {
+                               CPLFree(_src_srs);
+                               CPLFree(_dst_srs);
+                       }
 
                        return NULL;
                }
 
                switch (i) {
                        case 1:
-                               snprintf(
-                                       transform_opts[i],
-                                       sizeof(char) * (strlen("DST_SRS=") + strlen(_dst_srs) + 1),
-                                       "DST_SRS=%s",
-                                       _dst_srs
-                               );
+                               if (_dst_srs != NULL) {
+                                       snprintf(
+                                               transform_opts[i],
+                                               sizeof(char) * (strlen("DST_SRS=") + strlen(_dst_srs) + 1),
+                                               "DST_SRS=%s",
+                                               _dst_srs
+                                       );
+                               }
+                               else
+                                       sprintf(transform_opts[i], "%s", "DST_SRS=");
                                break;
                        case 0:
-                               snprintf(
-                                       transform_opts[i],
-                                       sizeof(char) * (strlen("SRC_SRS=") + strlen(_src_srs) + 1),
-                                       "SRC_SRS=%s",
-                                       _src_srs
-                               );
+                               if (_src_srs != NULL) {
+                                       snprintf(
+                                               transform_opts[i],
+                                               sizeof(char) * (strlen("SRC_SRS=") + strlen(_src_srs) + 1),
+                                               "SRC_SRS=%s",
+                                               _src_srs
+                                       );
+                               }
+                               else
+                                       sprintf(transform_opts[i], "%s", "SRC_SRS=");
                                break;
                }
                RASTER_DEBUGF(4, "transform_opts[%d] = %s", i, transform_opts[i]);
        }
        transform_opts[transform_opts_len] = NULL;
-       CPLFree(_src_srs);
-       CPLFree(_dst_srs);
+       if (_src_srs != NULL) CPLFree(_src_srs);
 
        /* transformation object for building dst dataset */
        transform_arg = GDALCreateGenImgProjTransformer2(src_ds, NULL, transform_opts);
@@ -9384,6 +9414,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9399,6 +9430,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9436,6 +9468,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9469,6 +9502,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9535,6 +9569,7 @@ rt_raster rt_raster_gdal_warp(
 
                        for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                        rtdealloc(transform_opts);
+                       if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                        return NULL;
                }
@@ -9565,6 +9600,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9605,6 +9641,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9627,6 +9664,7 @@ rt_raster rt_raster_gdal_warp(
 
                        for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                        rtdealloc(transform_opts);
+                       if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                        return NULL;
                }
@@ -9661,6 +9699,7 @@ rt_raster rt_raster_gdal_warp(
 
                                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                                rtdealloc(transform_opts);
+                               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                                return NULL;
                        }
@@ -9678,6 +9717,7 @@ rt_raster rt_raster_gdal_warp(
 
                                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                                rtdealloc(transform_opts);
+                               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                                return NULL;
                        }
@@ -9705,6 +9745,7 @@ rt_raster rt_raster_gdal_warp(
 
                                                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                                                rtdealloc(transform_opts);
+                                               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                                                return NULL;
                                        }
@@ -9734,6 +9775,7 @@ rt_raster rt_raster_gdal_warp(
 
                                                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                                                rtdealloc(transform_opts);
+                                               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                                                return NULL;
                                        }
@@ -9784,6 +9826,7 @@ rt_raster rt_raster_gdal_warp(
 
                                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                                rtdealloc(transform_opts);
+                               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                                return NULL;
                        }
@@ -9814,6 +9857,7 @@ rt_raster rt_raster_gdal_warp(
 
                                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                                rtdealloc(transform_opts);
+                               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                                return NULL;
                        }
@@ -9842,6 +9886,7 @@ rt_raster rt_raster_gdal_warp(
 
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
 
                return NULL;
        }
@@ -9865,6 +9910,38 @@ rt_raster rt_raster_gdal_warp(
 
                GDALClose(src_ds);
 
+               for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
+               rtdealloc(transform_opts);
+               if (_dst_srs != NULL) CPLFree(_dst_srs);
+
+               return NULL;
+       }
+
+       /* set dst srs */
+       if (_dst_srs != NULL) {
+               cplerr = GDALSetProjection(dst_ds, _dst_srs);
+               CPLFree(_dst_srs);
+               if (cplerr != CE_None) {
+                       rterror("rt_raster_gdal_warp: Unable to set projection");
+
+                       GDALClose(dst_ds);
+                       GDALClose(src_ds);
+
+                       for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
+                       rtdealloc(transform_opts);
+
+                       return NULL;
+               }
+       }
+
+       /* set dst geotransform */
+       cplerr = GDALSetGeoTransform(dst_ds, _gt);
+       if (cplerr != CE_None) {
+               rterror("rt_raster_gdal_warp: Unable to set geotransform");
+
+               GDALClose(dst_ds);
+               GDALClose(src_ds);
+
                for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
                rtdealloc(transform_opts);
 
@@ -9929,36 +10006,6 @@ rt_raster rt_raster_gdal_warp(
                }
        }
 
-       /* set dst srs */
-       _dst_srs = rt_util_gdal_convert_sr((NULL == dst_srs ? src_srs : dst_srs), 1);
-       cplerr = GDALSetProjection(dst_ds, _dst_srs);
-       CPLFree(_dst_srs);
-       if (cplerr != CE_None) {
-               rterror("rt_raster_gdal_warp: Unable to set projection");
-
-               GDALClose(dst_ds);
-               GDALClose(src_ds);
-
-               for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
-               rtdealloc(transform_opts);
-
-               return NULL;
-       }
-
-       /* set dst geotransform */
-       cplerr = GDALSetGeoTransform(dst_ds, _gt);
-       if (cplerr != CE_None) {
-               rterror("rt_raster_gdal_warp: Unable to set geotransform");
-
-               GDALClose(dst_ds);
-               GDALClose(src_ds);
-
-               for (i = 0; i < transform_opts_len; i++) rtdealloc(transform_opts[i]);
-               rtdealloc(transform_opts);
-
-               return NULL;
-       }
-
        /* create transformation object */
        transform_arg = imgproj_arg = GDALCreateGenImgProjTransformer2(src_ds, dst_ds, transform_opts);
        if (NULL == transform_arg) {
index 71e209e87d6d657faf49efea766816412e43b735..7c9fc4c486509e2cf9897f7f60da03b3e731f10b 100644 (file)
@@ -1594,8 +1594,9 @@ rt_raster rt_raster_from_gdal_dataset(GDALDatasetH ds);
  *
  * @return the warped raster or NULL
  */
-rt_raster rt_raster_gdal_warp(rt_raster raster, const char *src_srs,
-       const char *dst_srs,
+rt_raster rt_raster_gdal_warp(
+       rt_raster raster,
+       const char *src_srs, const char *dst_srs,
        double *scale_x, double *scale_y,
        int *width, int *height,
        double *ul_xw, double *ul_yw,
index 1c58ee43f87ec126c9045b5c3a5f3030c7472fed..f0da14e22ac25ebe12e63fed2a26590e3d50d525 100644 (file)
@@ -11259,17 +11259,21 @@ Datum RASTER_GDALWarp(PG_FUNCTION_ARGS)
                dst_srid = src_srid;
        POSTGIS_RT_DEBUGF(4, "destination SRID: %d", dst_srid);
 
-       /* source SRID = SRID_UNKNOWN */
-       if (src_srid == SRID_UNKNOWN) {
-               /* target SRID != src SRID, error */
-               if (dst_srid != src_srid) {
-                       elog(ERROR, "RASTER_GDALWarp: Input raster has unknown (%d) SRID", src_srid);
-                       rt_raster_destroy(raster);
-                       PG_FREE_IF_COPY(pgraster, 0);
-                       PG_RETURN_NULL();
+       /* target SRID != src SRID, error */
+       if (src_srid == SRID_UNKNOWN && dst_srid != src_srid) {
+               elog(ERROR, "RASTER_GDALWarp: Input raster has unknown (%d) SRID", src_srid);
+               rt_raster_destroy(raster);
+               PG_FREE_IF_COPY(pgraster, 0);
+               PG_RETURN_NULL();
+       }
+       /* target SRID == src SRID, no reprojection */
+       else if (dst_srid == src_srid) {
+               /* set geotransform */
+               if (src_srid == SRID_UNKNOWN) {
+                       double gt[6] = {0, 10, 0, 0, 0, -10};
+                       rt_raster_set_geotransform_matrix(raster, gt);
                }
 
-               /* target SRID == src SRID, special */
                no_srid = 1;
        }
 
@@ -11368,12 +11372,8 @@ Datum RASTER_GDALWarp(PG_FUNCTION_ARGS)
        }
 
        /* get srses from srids */
-       /* source srs */
-       /* no SRID, use EPSG:4326 (WGS84) */
-       if (no_srid)
-               src_srs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs";
-       /* legitimate SRID */
-       else {
+       if (!no_srid) {
+               /* source srs */
                src_srs = rtpg_getSR(src_srid);
                if (NULL == src_srs) {
                        elog(ERROR, "RASTER_GDALWarp: Input raster has unknown SRID (%d)", src_srid);
@@ -11381,19 +11381,12 @@ Datum RASTER_GDALWarp(PG_FUNCTION_ARGS)
                        PG_FREE_IF_COPY(pgraster, 0);
                        PG_RETURN_NULL();
                }
-       }
-       POSTGIS_RT_DEBUGF(4, "src srs: %s", src_srs);
+               POSTGIS_RT_DEBUGF(4, "src srs: %s", src_srs);
 
-       /* target srs */
-       /* no SRID, use src_srs */
-       if (no_srid)
-               dst_srs = src_srs;
-       /* legitimate SRID */
-       else if (dst_srid != SRID_UNKNOWN) {
                dst_srs = rtpg_getSR(dst_srid);
                if (NULL == dst_srs) {
                        elog(ERROR, "RASTER_GDALWarp: Target SRID (%d) is unknown", dst_srid);
-                       if (!no_srid && NULL != src_srs) pfree(src_srs);
+                       if (!no_srid) pfree(src_srs);
                        rt_raster_destroy(raster);
                        PG_FREE_IF_COPY(pgraster, 0);
                        PG_RETURN_NULL();
@@ -11413,17 +11406,21 @@ Datum RASTER_GDALWarp(PG_FUNCTION_ARGS)
        rt_raster_destroy(raster);
        PG_FREE_IF_COPY(pgraster, 0);
        if (!no_srid) {
-               if (NULL != src_srs) pfree(src_srs);
-               if (NULL != dst_srs) pfree(dst_srs);
+               pfree(src_srs);
+               pfree(dst_srs);
        }
        if (!rast) {
                elog(ERROR, "RASTER_band: Could not create transformed raster");
                PG_RETURN_NULL();
        }
 
-       /* add target SRID but only if we're not using a no SRID */
-       if (!no_srid)
-               rt_raster_set_srid(rast, dst_srid);
+       /* add target SRID */
+       rt_raster_set_srid(rast, dst_srid);
+
+       if (no_srid && src_srid == SRID_UNKNOWN) {
+               double gt[6] = {0, 1, 0, 0, 0, -1};
+               rt_raster_set_geotransform_matrix(rast, gt);
+       }
 
        pgrast = rt_raster_serialize(rast);
        rt_raster_destroy(rast);
index fb05c02447eaf772a9f5c9abb67cb9e5e5933978..76593ee06efb8233809b3ef3e4a53b367e792f26 100644 (file)
@@ -247,7 +247,30 @@ INSERT INTO raster_gdalwarp_dst (rid, rast) VALUES (
                ST_SetSRID(rast, 0),
                'NearestNeighbor', 0.125,
                NULL,
-               500., 500.
+               NULL, NULL,
+               NULL, NULL,
+               NULL, NULL,
+               5, 5
+       ) FROM raster_gdalwarp_src)
+), (
+       0.26, (SELECT _st_gdalwarp(
+               ST_SetSRID(rast, 0),
+               'NearestNeighbor', 0.125,
+               NULL,
+               NULL, NULL,
+               NULL, NULL,
+               NULL, NULL,
+               2, 2
+       ) FROM raster_gdalwarp_src)
+), (
+       0.27, (SELECT _st_gdalwarp(
+               ST_SetSRID(rast, 0),
+               'NearestNeighbor', 0.125,
+               NULL,
+               NULL, NULL,
+               NULL, NULL,
+               NULL, NULL,
+               100, 100
        ) FROM raster_gdalwarp_src)
 );
 
index ca376cc32251739ee1262a4d811c43c5ab4e3df4..9ce521b093a1fc59191c496efd069c681e10da63 100644 (file)
@@ -22,7 +22,9 @@ NOTICE:  Values must be provided for both X and Y when specifying the scale.  Re
 0.22|993310|26|26|1|500.000|500.000|0.000|6.000|950452.000|1396632.000|t|t|t
 0.23|984269|12|8|1|0.012|-0.012|0.000|0.000|-107.029|50.206|t|t|t
 0.24|974269|12|8|1|0.012|-0.012|0.000|0.000|-107.029|50.206|t|t|t
-0.25|0|20|20|1|500.000|500.000|0.000|0.000|-500000.000|590000.000|t|t|t
+0.25|0|5|5|1|1.000|-1.000|0.000|0.000|0.000|0.000|t|t|t
+0.26|0|2|2|1|1.000|-1.000|0.000|0.000|0.000|0.000|t|t|t
+0.27|0|100|100|1|1.000|-1.000|0.000|0.000|0.000|0.000|t|t|t
 0.3|994269|12|8|1|0.012|-0.012|0.000|0.000|-107.029|50.206|t|t|t
 0.4|993310|24|24|1|500.000|500.000|0.000|0.000|950732.188|1397281.783|t|t|t
 0.5|992163|10|10|1|1000.000|-1000.000|0.000|0.000|-500000.000|600000.000|t|t|t