How to use the gdal.AutoCreateWarpedVRT function in GDAL

To help you get started, we’ve selected a few GDAL examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github giohappy / gdal2cesium / gdal2cesium.py View on Github external
"Either gdal2tiles with parameter -p 'raster' or use another GIS software for georeference e.g. gdal_transform -gcp / -a_ullr / -a_srs")
        
        in_srs_code = self.in_srs.GetAttrValue("AUTHORITY", 0)
        in_ds_srs = osr.SpatialReference()
        res = in_ds_srs.ImportFromWkt(in_ds.GetProjection())
        
        if res != 0 and in_srs_code is None:
                print "ERROR! The inumpyut file %s has no SRS associated and no SRS has been defined in inumpyut (-s parameter)" % _inumpyut
                exit(1)
 
        if self.in_srs:
            if in_ds_srs.ExportToProj4() != self.out_srs.ExportToProj4():
                if (self.in_srs.ExportToProj4() != self.out_srs.ExportToProj4()) or (in_ds.GetGCPCount() != 0):
                    print "WARNING! Inumpyut file %s has a SR different from EPSG:4326 (WGS84). This can make the processing significantly slow." % _inumpyut
                    # Generation of VRT dataset in tile projection, default 'nearest neighbour' warping
                    out_ds = gdal.AutoCreateWarpedVRT( in_ds, self.in_srs_wkt, self.out_srs.ExportToWkt() )

                    # TODO: HIGH PRIORITY: Correction of AutoCreateWarpedVRT according the max zoomlevel for correct direct warping!!!

                    if self.options.verbose:
                        print("Warping of the raster by AutoCreateWarpedVRT (result saved into 'tiles.vrt')")
                    out_ds.GetDriver().CreateCopy("%s.vrt" % _inumpyut, out_ds)
                    inumpyut_or_vrt = "%s.vrt" % _inumpyut

                    # Note: self.in_srs and self.in_srs_wkt contain still the non-warped reference system!!!

        else:
            self.error("Inumpyut file has unknown SRS.", "Use --s_srs ESPG:xyz (or similar) to provide source reference system." )

        if out_ds and self.options.verbose:
            print("Projected file:", "tiles.vrt", "( %sP x %sL - %s bands)" % (out_ds.RasterXSize, out_ds.RasterYSize, out_ds.RasterCount))
github buddebej / dem2tiles / tiler-tools / tiler_backend.py View on Github external
def calc_zoom(self, zoom_parm):
        'determine and set a list of zoom levels to generate'
    #----------------------------

        # check raster parameters to find default zoom range
        ld('automatic zoom levels')

        # modify target srs to allow charts crossing meridian 180
        shifted_srs = self.shift_srs()

        t_ds = gdal.AutoCreateWarpedVRT(self.src_ds, None, txt2wkt(shifted_srs))
        geotr = t_ds.GetGeoTransform()
        res = (geotr[1], geotr[5])
        max_zoom = max(self.res2zoom_xy(res))

        # calculate min_zoom
        ul_c = (geotr[0], geotr[3])
        lr_c = gdal.ApplyGeoTransform(geotr, t_ds.RasterXSize, t_ds.RasterYSize)
        wh = (lr_c[0]-ul_c[0], ul_c[1]-lr_c[1])
        ld('ul_c, lr_c, wh', ul_c, lr_c, wh)
        min_zoom = min(self.res2zoom_xy([wh[i]/abs(self.tile_dim[i]) for i in (0, 1)]))

        self.set_zoom_range(zoom_parm, (min_zoom, max_zoom))
github nansencenter / nansat / nansat / vrt.py View on Github external
"""
        # VRT to be warped
        src_vrt = self.copy()

        # if destination GCPs are given: create and add fake GCPs to src
        dst_wkt = src_vrt._set_fake_gcps(dst_srs, dst_gcps, skip_gcps)

        if resize_only:
            src_vrt._set_geotransform_for_resize()
        else:
            src_vrt._set_gcps_geolocation_geotransform()

        # create Warped VRT GDAL Dataset
        warped_dataset = gdal.AutoCreateWarpedVRT(src_vrt.dataset, None, dst_wkt, resample_alg)

        # check if Warped VRT was created
        if warped_dataset is None:
            raise NansatGDALError('Cannot create warpedVRT!')

        # create VRT object from Warped VRT GDAL Dataset
        warped_vrt = VRT.copy_dataset(warped_dataset)

        # set x/y size, geo_transform, block_size
        warped_vrt._update_warped_vrt_xml(x_size, y_size, geo_transform, block_size, working_data_type)

        # apply thin-spline-transformation option
        if self.tps:
            warped_vrt.write_xml(warped_vrt.xml.replace('GCPTransformer', 'TPSTransformer'))

        # if given, add dst GCPs
github buddebej / ol3-dem / dem-preprocessing / tiler-tools / tiler_backend.py View on Github external
def calc_zoom(self, zoom_parm):
        'determine and set a list of zoom levels to generate'
    #----------------------------

        # check raster parameters to find default zoom range
        ld('automatic zoom levels')

        # modify target srs to allow charts crossing meridian 180
        shifted_srs = self.shift_srs()

        t_ds = gdal.AutoCreateWarpedVRT(self.src_ds, None, txt2wkt(shifted_srs))
        geotr = t_ds.GetGeoTransform()
        res = (geotr[1], geotr[5])
        max_zoom = max(self.res2zoom_xy(res))

        # calculate min_zoom
        ul_c = (geotr[0], geotr[3])
        lr_c = gdal.ApplyGeoTransform(geotr, t_ds.RasterXSize, t_ds.RasterYSize)
        wh = (lr_c[0]-ul_c[0], ul_c[1]-lr_c[1])
        ld('ul_c, lr_c, wh', ul_c, lr_c, wh)
        min_zoom = min(self.res2zoom_xy([wh[i]/abs(self.tile_dim[i]) for i in (0, 1)]))

        self.set_zoom_range(zoom_parm, (min_zoom, max_zoom))
github commenthol / gdal2tiles-leaflet / gdal2tiles.py View on Github external
1.0,
                ) and self.in_ds.GetGCPCount() == 0:
                self.error("There is no georeference - neither affine transformation (worldfile) nor GCPs. You can generate only 'raster' profile tiles."
                           ,
                           "Either gdal2tiles with parameter -p 'raster' or use another GIS software for georeference e.g. gdal_transform -gcp / -a_ullr / -a_srs"
                           )

            if self.in_srs:

                if self.in_srs.ExportToProj4() \
                    != self.out_srs.ExportToProj4() \
                    or self.in_ds.GetGCPCount() != 0:

                    # Generation of VRT dataset in tile projection, default 'nearest neighbour' warping

                    self.out_ds = gdal.AutoCreateWarpedVRT(self.in_ds,
                            self.in_srs_wkt, self.out_srs.ExportToWkt())

                    # TODO: HIGH PRIORITY: Correction of AutoCreateWarpedVRT according the max zoomlevel for correct direct warping!!!

                    if self.options.verbose:
                        print "Warping of the raster by AutoCreateWarpedVRT (result saved into 'tiles.vrt')"
                        self.out_ds.GetDriver().CreateCopy('tiles.vrt',
                                self.out_ds)

                    # Note: self.in_srs and self.in_srs_wkt contain still the non-warped reference system!!!

                    # Correction of AutoCreateWarpedVRT for NODATA values

                    if self.in_nodata != []:
                        import tempfile
                        tempfilename = tempfile.mktemp('-gdal2tiles.vrt'
github buddebej / dem2tiles / tiler-tools / tiler_backend.py View on Github external
# calculate zoom range
        self.calc_zoom(zoom_parm)
        self.max_zoom = self.zoom_range[0]

        # shift target SRS to avoid crossing 180 meridian
        shifted_srs = self.shift_srs(self.max_zoom)
        shift_x = GdalTransformer(SRC_SRS=shifted_srs, DST_SRS=self.proj_srs).transform_point((0, 0))[0]
        if shift_x != 0:
            self.proj_srs = shifted_srs
            self.proj2geog = GdalTransformer(SRC_SRS=self.proj_srs, DST_SRS=self.geog_srs)
            self.pix_origin = (self.pix_origin[0]-shift_x, self.pix_origin[1])
            self.tile_origin = (self.tile_origin[0]-shift_x, self.tile_origin[1])
            ld('new_srs', shifted_srs, 'shift_x', shift_x, 'pix_origin', self.pix_origin)

        # get corners at the target SRS
        target_ds = gdal.AutoCreateWarpedVRT(self.src_ds, None, txt2wkt(shifted_srs))
        target_bounds = GdalTransformer(target_ds).transform([
            (0, 0),
            (target_ds.RasterXSize, target_ds.RasterYSize)])

        # self.bounds are set to a world raster, now clip to the max tileset area
        self.bounds = ((target_bounds[0][0],
                      min(self.bounds[0][1], target_bounds[0][1])),
                     (target_bounds[1][0],
                      max(self.bounds[1][1], target_bounds[1][1])))

        ld('target raster')
        ld('Upper Left', self.bounds[0], target_bounds[0], self.proj2geog.transform([self.bounds[0], target_bounds[0]]))
        ld('Lower Right', self.bounds[1], target_bounds[1], self.proj2geog.transform([self.bounds[1], target_bounds[1]]))
#        orig_ul = GdalTransformer(SRC_SRS=self.geog_srs, DST_SRS=self.srs).transform_point(
#            self.proj2geog.transform_point(target_bounds[0]))
#        ld(orig_ul[0]-target_bounds[0][0], orig_ul)
github GitHubRGI / geopackage-python / Tiling / gdal2tiles_parallel.py View on Github external
if (self.in_ds.GetGeoTransform() ==
                (0.0, 1.0, 0.0, 0.0, 0.0, 1.0)) and (
                    self.in_ds.GetGCPCount() == 0):
                self.error(
                    "There is no georeference - neither affine transformation (worldfile) nor GCPs. You can generate only 'raster' profile tiles.",
                    "Either gdal2tiles with parameter -p 'raster' or use another GIS software for georeference e.g. gdal_transform -gcp / -a_ullr / -a_srs")

            if self.in_srs:

                if (self.in_srs.ExportToProj4() !=
                        self.out_srs.ExportToProj4()) or (
                            self.in_ds.GetGCPCount() != 0):

                    # Generation of VRT dataset in tile projection, default 'nearest neighbour' warping
                    self.out_ds = gdal.AutoCreateWarpedVRT(
                        self.in_ds, self.in_srs_wkt,
                        self.out_srs.ExportToWkt())

                    # TODO: HIGH PRIORITY: Correction of AutoCreateWarpedVRT according the max zoomlevel for correct direct warping!!!

                    if self.options.verbose:
                        print(
                            "Warping of the raster by AutoCreateWarpedVRT (result saved into 'tiles.vrt')")
                        self.out_ds.GetDriver().CreateCopy("tiles.vrt",
                                                           self.out_ds)

                        # Note: self.in_srs and self.in_srs_wkt contain still the non-warped reference system!!!

                        # Correction of AutoCreateWarpedVRT for NODATA values
                    if self.in_nodata != []:
                        tempfilename = mktemp('-gdal2tiles.vrt')
github OpenDataAnalytics / gaia / gaia / geo / gdal_functions.py View on Github external
"""
    # Open source dataset
    src_ds = get_dataset(src)

    # Define target SRS
    dst_srs = osr.SpatialReference()
    dst_srs.ImportFromEPSG(int(epsg))
    dst_wkt = dst_srs.ExportToWkt()

    # Resampling might be passed as a string
    if not isinstance(resampling, int):
        resampling = getattr(gdal, resampling)

    # Call AutoCreateWarpedVRT() to fetch default values
    # for target raster dimensions and geotransform
    reprojected_ds = gdal.AutoCreateWarpedVRT(src_ds,
                                              None,
                                              dst_wkt,
                                              resampling,
                                              error_threshold)

    # Create the final warped raster
    if dst:
        gdal.GetDriverByName('GTiff').CreateCopy(dst, reprojected_ds)
    return reprojected_ds