Python gdal_array.NumericTypeCodeToGDALTypeCode方法代碼示例

# 需要導入模塊: from osgeo import gdal_array [as 別名]
# 或者: from osgeo.gdal_array import NumericTypeCodeToGDALTypeCode [as 別名]
def to_tiff(self, filename):
        '''
        geo.to_tiff(filename)

        Saves GeoRaster as geotiff filename.tif with type datatype

        If GeoRaster does not have datatype, then it tries to assign a type.
        You can assign the type yourself by setting
         geo.datatype = 'gdal.GDT_'+type
        '''
        if self.datatype is None:
            self.datatype = gdal_array.NumericTypeCodeToGDALTypeCode(self.raster.data.dtype)
            if self.datatype is None:
                if self.raster.data.dtype.name.find('int') !=- 1:
                    self.raster = self.raster.astype(np.int32)
                    self.datatype = gdal_array.NumericTypeCodeToGDALTypeCode(self.raster.data.dtype)
                else:
                    self.raster = self.raster.astype(np.float64)
                    self.datatype = gdal_array.NumericTypeCodeToGDALTypeCode(self.raster.data.dtype)
        self.raster.data[self.raster.mask] = self.nodata_value
        create_geotiff(filename, self.raster, gdal.GetDriverByName('GTiff'), self.nodata_value,
                       self.shape[1], self.shape[0], self.geot, self.projection, self.datatype) 

# 需要導入模塊: from osgeo import gdal_array [as 別名]
# 或者: from osgeo.gdal_array import NumericTypeCodeToGDALTypeCode [as 別名]
def save_array_as_image(array, path, geotransform, projection, format = "GTiff"):
    """
    Saves a given array as a geospatial image to disk in the format 'format'. The datatype will be of one corresponding
    to
    Array must be gdal format: [bands, y, x].

    Parameters
    ----------
    array
        A Numpy array containing the values to be saved to a raster
    path
        The path to the location to save the output raster to
    geotransform
        The geotransform of the image to be saved. See note.
    projection
        The projection, as wkt, of the image to be saved. See note.
    format

    Returns
    -------

    """
    driver = gdal.GetDriverByName(format)
    type_code = gdal_array.NumericTypeCodeToGDALTypeCode(array.dtype)
    # If array is 2d, give it an extra dimension.
    if len(array.shape) == 2:
        array = np.expand_dims(array, axis=0)
    out_dataset = driver.Create(
        path,
        xsize=array.shape[2],
        ysize=array.shape[1],
        bands=array.shape[0],
        eType=type_code
    )
    out_dataset.SetGeoTransform(geotransform)
    out_dataset.SetProjection(projection)
    out_array = out_dataset.GetVirtualMemArray(eAccess=gdal.GA_Update).squeeze()
    out_array[...] = array
    out_array = None
    out_dataset = None
    return path 

# 需要導入模塊: from osgeo import gdal_array [as 別名]
# 或者: from osgeo.gdal_array import NumericTypeCodeToGDALTypeCode [as 別名]
def writeGTiff(dstImage, outputPath, projection, coordinates, dtype = None):
    """
    @brief Writes a GeoTiff file created from an existing coordinate-system
    @Note Coordinates: [Top-Left X, W-E Resolution, 0, Top Left Y, 0, N-S Resolution]
    """
    driver = gdal.GetDriverByName('GTiff')
    #Add dimension for a single band-image
    if(dstImage.ndim is 2):
        dstImage = dstImage[..., np.newaxis]  
    #Set output dtype if not specified
    if(dtype == None):
        dtype = gdal_array.NumericTypeCodeToGDALTypeCode(dstImage.dtype)
    dataset = driver.Create(outputPath, dstImage.shape[1], dstImage.shape[0] , dstImage.shape[2], dtype)
    if(dataset == None):
        raise OSError ("GDAL Could not create {0}".format(outputPath))
    if(not coordinates or len(coordinates) != 6):
        raise ValueError("Geotransform empty or unknown: {0}".format(coordinates))
    dataset.SetGeoTransform(coordinates)  
    if(not projection):
        raise ValueError("Projection empty or unknown: {0}".format(projection))
    dataset.SetProjection(projection)
    for index in range(dstImage.shape[2]):
        dataset.GetRasterBand(index+1).WriteArray(dstImage[:,:,index])
    dataset.FlushCache()  # Write to disk.
    dataset = None
    driver = None
    return 0 

# 需要導入模塊: from osgeo import gdal_array [as 別名]
# 或者: from osgeo.gdal_array import NumericTypeCodeToGDALTypeCode [as 別名]
def write_output(raster, output, image_ds, gdal_frmt, ndv, band_names=None):
    """ Write raster to output file """
    from osgeo import gdal, gdal_array

    logger.debug('Writing output to disk')

    driver = gdal.GetDriverByName(str(gdal_frmt))

    if len(raster.shape) > 2:
        nband = raster.shape[2]
    else:
        nband = 1

    ds = driver.Create(
        output,
        image_ds.RasterXSize, image_ds.RasterYSize, nband,
        gdal_array.NumericTypeCodeToGDALTypeCode(raster.dtype.type)
    )

    if band_names is not None:
        if len(band_names) != nband:
            logger.error('Did not get enough names for all bands')
            sys.exit(1)

    if raster.ndim > 2:
        for b in range(nband):
            logger.debug('    writing band {b}'.format(b=b + 1))
            ds.GetRasterBand(b + 1).WriteArray(raster[:, :, b])
            ds.GetRasterBand(b + 1).SetNoDataValue(ndv)

            if band_names is not None:
                ds.GetRasterBand(b + 1).SetDescription(band_names[b])
                ds.GetRasterBand(b + 1).SetMetadata({
                    'band_{i}'.format(i=b + 1): band_names[b]
                })
    else:
        logger.debug('    writing band')
        ds.GetRasterBand(1).WriteArray(raster)
        ds.GetRasterBand(1).SetNoDataValue(ndv)

        if band_names is not None:
            ds.GetRasterBand(1).SetDescription(band_names[0])
            ds.GetRasterBand(1).SetMetadata({'band_1': band_names[0]})

    ds.SetProjection(image_ds.GetProjection())
    ds.SetGeoTransform(image_ds.GetGeoTransform())

    ds = None


# RESULT UTILITIES 

# 需要導入模塊: from osgeo import gdal_array [as 別名]
# 或者: from osgeo.gdal_array import NumericTypeCodeToGDALTypeCode [as 別名]
def create_raster_dataset(data, coords, projection=None, nodata=-9999):
    """ Create In-Memory Raster Dataset

    Parameters
    ----------
    data : :class:`numpy:numpy.ndarray`
        Array of shape (rows, cols) or (bands, rows, cols) containing
        the data values.
    coords : :class:`numpy:numpy.ndarray`
        Array of shape (nrows, ncols, 2) containing pixel center coordinates
        or
        Array of shape (nrows+1, ncols+1, 2) containing pixel edge coordinates
    projection : osr object
        Spatial reference system of the used coordinates, defaults to None.
    nodata : int
        Value of NODATA

    Returns
    -------
    dataset : gdal.Dataset
        In-Memory raster dataset

    Note
    ----
    The origin of the provided data and coordinates is UPPER LEFT.
    """

    # align data
    data = data.copy()
    if data.ndim == 2:
        data = data[np.newaxis, ...]
    bands, rows, cols = data.shape

    # create In-Memory Raster with correct dtype
    mem_drv = gdal.GetDriverByName("MEM")
    gdal_type = gdal_array.NumericTypeCodeToGDALTypeCode(data.dtype)
    dataset = mem_drv.Create("", cols, rows, bands, gdal_type)

    # initialize geotransform
    x_ps, y_ps = coords[1, 1] - coords[0, 0]
    if data.shape[-2:] == coords.shape[0:2]:
        upper_corner_x = coords[0, 0, 0] - x_ps / 2
        upper_corner_y = coords[0, 0, 1] - y_ps / 2
    else:
        upper_corner_x = coords[0, 0, 0]
        upper_corner_y = coords[0, 0, 1]
    geotran = [upper_corner_x, x_ps, 0, upper_corner_y, 0, y_ps]
    dataset.SetGeoTransform(geotran)

    if projection:
        dataset.SetProjection(projection.ExportToWkt())

    # set np.nan to nodata
    dataset.GetRasterBand(1).SetNoDataValue(nodata)

    for i, band in enumerate(data, start=1):
        dataset.GetRasterBand(i).WriteArray(band)

    return dataset