build_overview()

static int build_overview ( int  idx, RTLOADERCFG *  config, RASTERINFO *  info, uint32_t  ovx, STRINGBUFFER *  tileset, STRINGBUFFER *  buffer  )

static

Definition at line 1340 of file raster2pgsql.c.

                                                                                                                          {
        GDALDatasetH hdsSrc;
        VRTDatasetH hdsOv;
        VRTSourcedRasterBandH hbandOv;
        double gtOv[6] = {0.};
        int dimOv[2] = {0};
 
        uint32_t j = 0;
        int factor;
        const char *ovtable = NULL;
 
        VRTDatasetH hdsDst;
        VRTSourcedRasterBandH hbandDst;
        int tile_size[2] = {0};
        int _tile_size[2] = {0};
        int ntiles[2] = {1, 1};
        int xtile = 0;
        int ytile = 0;
        double gt[6] = {0.};
 
        rt_raster rast = NULL;
        char *hex;
        uint32_t hexlen = 0;
 
        hdsSrc = GDALOpenShared(config->rt_file[idx], GA_ReadOnly);
        if (hdsSrc == NULL) {
                rterror(_("build_overview: Could not open raster: %s"), config->rt_file[idx]);
                return 0;
        }
 
        /* working copy of geotransform matrix */
        memcpy(gtOv, info->gt, sizeof(double) * 6);
 
        if (ovx >= config->overview_count) {
                rterror(_("build_overview: Invalid overview index: %d"), ovx);
                return 0;
        }
        factor = config->overview[ovx];
        ovtable = (const char *) config->overview_table[ovx];
 
        /* factor must be within valid range */
        if (factor < MINOVFACTOR || factor > MAXOVFACTOR) {
                rterror(_("build_overview: Overview factor %d is not between %d and %d"), factor, MINOVFACTOR, MAXOVFACTOR);
                return 0;
        }
 
        dimOv[0] = (int) (info->dim[0] + (factor / 2)) / factor;
        dimOv[1] = (int) (info->dim[1] + (factor / 2)) / factor;
 
        /* create VRT dataset */
        hdsOv = VRTCreate(dimOv[0], dimOv[1]);
        /*
        GDALSetDescription(hdsOv, "/tmp/ov.vrt");
        */
        GDALSetProjection(hdsOv, info->srs);
 
        /* adjust scale */
        gtOv[1] *= factor;
        gtOv[5] *= factor;
 
        GDALSetGeoTransform(hdsOv, gtOv);
 
        /* add bands as simple sources */
        for (j = 0; j < info->nband_count; j++) {
                GDALAddBand(hdsOv, info->gdalbandtype[j], NULL);
                hbandOv = (VRTSourcedRasterBandH) GDALGetRasterBand(hdsOv, j + 1);
 
                if (info->hasnodata[j])
                        GDALSetRasterNoDataValue(hbandOv, info->nodataval[j]);
 
                VRTAddSimpleSource(
                        hbandOv, GDALGetRasterBand(hdsSrc, info->nband[j]),
                        0, 0,
                        info->dim[0], info->dim[1],
                        0, 0,
                        dimOv[0], dimOv[1],
                        "near", VRT_NODATA_UNSET
                );
        }
 
        /* make sure VRT reflects all changes */
        VRTFlushCache(hdsOv);
 
        /* decide on tile size */
        if (!config->tile_size[0])
                tile_size[0] = dimOv[0];
        else
                tile_size[0] = config->tile_size[0];
        if (!config->tile_size[1])
                tile_size[1] = dimOv[1];
        else
                tile_size[1] = config->tile_size[1];
 
        /* number of tiles */
        if (
                tile_size[0] != dimOv[0] &&
                tile_size[1] != dimOv[1]
        ) {
                ntiles[0] = (dimOv[0] + tile_size[0] -  1) / tile_size[0];
                ntiles[1] = (dimOv[1] + tile_size[1]  - 1) / tile_size[1];
        }
 
        /* working copy of geotransform matrix */
        memcpy(gt, gtOv, sizeof(double) * 6);
 
        /* tile overview */
        /* each tile is a VRT with constraints set for just the data required for the tile */
        for (ytile = 0; ytile < ntiles[1]; ytile++) {
 
                /* edge y tile */
                if (!config->pad_tile && ntiles[1] > 1 && (ytile + 1) == ntiles[1])
                        _tile_size[1] = dimOv[1] - (ytile * tile_size[1]);
                else
                        _tile_size[1] = tile_size[1];
 
                for (xtile = 0; xtile < ntiles[0]; xtile++) {
                        /*
                        char fn[100];
                        sprintf(fn, "/tmp/ovtile%d.vrt", (ytile * ntiles[0]) + xtile);
                        */
 
                        /* edge x tile */
                        if (!config->pad_tile && ntiles[0] > 1 && (xtile + 1) == ntiles[0])
                                _tile_size[0] = dimOv[0] - (xtile * tile_size[0]);
                        else
                                _tile_size[0] = tile_size[0];
 
                        /* compute tile's upper-left corner */
                        GDALApplyGeoTransform(
                                gtOv,
                                xtile * tile_size[0], ytile * tile_size[1],
                                &(gt[0]), &(gt[3])
                        );
 
                        /* create VRT dataset */
                        hdsDst = VRTCreate(_tile_size[0], _tile_size[1]);
                        /*
                        GDALSetDescription(hdsDst, fn);
                        */
                        GDALSetProjection(hdsDst, info->srs);
                        GDALSetGeoTransform(hdsDst, gt);
 
                        /* add bands as simple sources */
                        for (j = 0; j < info->nband_count; j++) {
                                GDALAddBand(hdsDst, info->gdalbandtype[j], NULL);
                                hbandDst = (VRTSourcedRasterBandH) GDALGetRasterBand(hdsDst, j + 1);
 
                                if (info->hasnodata[j])
                                        GDALSetRasterNoDataValue(hbandDst, info->nodataval[j]);
 
                                VRTAddSimpleSource(
                                        hbandDst, GDALGetRasterBand(hdsOv, j + 1),
                                        xtile * tile_size[0], ytile * tile_size[1],
                                        _tile_size[0], _tile_size[1],
                                        0, 0,
                                        _tile_size[0], _tile_size[1],
                                        "near", VRT_NODATA_UNSET
                                );
                        }
 
                        /* make sure VRT reflects all changes */
                        VRTFlushCache(hdsDst);
 
                        /* convert VRT dataset to rt_raster */
                        rast = rt_raster_from_gdal_dataset(hdsDst);
                        if (rast == NULL) {
                                rterror(_("build_overview: Could not convert VRT dataset to PostGIS raster"));
                                GDALClose(hdsDst);
                                return 0;
                        }
 
                        /* set srid if provided */
                        rt_raster_set_srid(rast, info->srid);
 
                        /* convert rt_raster to hexwkb */
                        hex = rt_raster_to_hexwkb(rast, FALSE, &hexlen);
                        raster_destroy(rast);
 
                        if (hex == NULL) {
                                rterror(_("build_overview: Could not convert PostGIS raster to hex WKB"));
                                GDALClose(hdsDst);
                                return 0;
                        }
 
                        /* add hexwkb to tileset */
                        append_stringbuffer(tileset, hex);
 
                        GDALClose(hdsDst);
 
                        /* flush if tileset gets too big */
                        if (tileset->length >= config->max_tiles_per_copy) {
                                if (!insert_records(
                                        config->schema, ovtable, config->raster_column,
                                        (config->file_column ? config->rt_filename[idx] : NULL), config->file_column_name,
                                        config->copy_statements, config->out_srid,
                                        tileset, buffer
                                )) {
                                        rterror(_("build_overview: Could not convert raster tiles into INSERT or COPY statements"));
                                        GDALClose(hdsSrc);
                                        return 0;
                                }
 
                                rtdealloc_stringbuffer(tileset, 0);
                        }
                }
        }
 
        GDALClose(hdsOv);
        GDALClose(hdsSrc);
        return 1;
}

References _, append_stringbuffer(), buffer(), raster_loader_config::copy_statements, rasterinfo_t::dim, FALSE, raster_loader_config::file_column, raster_loader_config::file_column_name, rasterinfo_t::gdalbandtype, rasterinfo_t::gt, rasterinfo_t::hasnodata, insert_records(), stringbuffer_t::length, raster_loader_config::max_tiles_per_copy, MAXOVFACTOR, MINOVFACTOR, rasterinfo_t::nband, rasterinfo_t::nband_count, rasterinfo_t::nodataval, raster_loader_config::out_srid, raster_loader_config::overview, raster_loader_config::overview_count, raster_loader_config::overview_table, raster_loader_config::pad_tile, raster_loader_config::raster_column, raster_destroy(), raster_loader_config::rt_file, raster_loader_config::rt_filename, rt_raster_from_gdal_dataset(), rt_raster_set_srid(), rt_raster_to_hexwkb(), rtdealloc_stringbuffer(), rterror(), raster_loader_config::schema, rasterinfo_t::srid, rasterinfo_t::srs, and raster_loader_config::tile_size.

Referenced by process_rasters().

Here is the call graph for this function:

Here is the caller graph for this function: