We should consider cell height and width when clipping rasters #1301

Closed
akbargumbira opened this Issue Oct 30, 2014 · 2 comments

Comments

Projects
None yet
5 participants
@akbargumbira
Member

akbargumbira commented Oct 30, 2014

Problem

When we do the clipping in the case when the hazard and the exposure are raster, we only consider the width of the cell size. See get_wgs84_resolution:

def get_wgs84_resolution(layer):
    """Return resolution of raster layer in EPSG:4326.

    If input layer is already in EPSG:4326, simply return the resolution
    If not, work it out based on EPSG:4326 representations of its extent.

    :param layer: Raster layer
    :type layer: QgsRasterLayer or QgsMapLayer

    :returns: The resolution of the given layer.
    :rtype: float

    """

    msg = tr(
        'Input layer to getWGS84resolution must be a raster layer. '
        'I got: %s' % str(layer.type())[1:-1])
    if not layer.type() == QgsMapLayer.RasterLayer:
        raise RuntimeError(msg)

    if layer.crs().authid() == 'EPSG:4326':
        cell_size = layer.rasterUnitsPerPixelX()
    else:
        # Otherwise, work it out based on EPSG:4326 representations of
        # its extent

        # Reproject extent to EPSG:4326
        geo_crs = QgsCoordinateReferenceSystem()
        geo_crs.createFromSrid(4326)
        transform = QgsCoordinateTransform(layer.crs(), geo_crs)
        extent = layer.extent()
        projected_extent = transform.transformBoundingBox(extent)

        # Estimate cell size
        columns = layer.width()
        geo_width = abs(
            projected_extent.xMaximum() -
            projected_extent.xMinimum())
        cell_size = geo_width / columns

    return cell_size

When we do the clipping, to get the resolution of the clipped layer, here is our logic:

if hazard_geo_cell_size < exposure_geo_cell_size and \
                        allow_resampling_flag:
                    cell_size = hazard_geo_cell_size
                    layer_extent = hazard_geoextent
                else:
                    cell_size = exposure_geo_cell_size
                    layer_extent = exposure_geoextent

It means that we will always get square cell of the clipped layer. If both input layers are rectangular, e.g have these resolutions: Hazard (7, 2), Exposure (10, 1), the clipped layer will have the cell size (7, 7).

This also will make other function not working as we expect e.g for the case when hazard is raster and the exposure is vector, we will buffer the extent to make the interpolation work for edge pixels using:

 adjusted_geo_extent = get_buffered_extent(
                    geo_extent,
                    hazard_geo_cell_size)

which the buffered_extent itself is expecting the cell size as a tuple (see function buffered_bounding_box in safe.storage.utilities:

 bbox = copy.copy(list(bbox))

    if resolution is None:
        return bbox

    try:
        resx, resy = resolution
    except TypeError:
        resx = resy = resolution

    bbox[0] -= resx
    bbox[1] -= resy
    bbox[2] += resx
    bbox[3] += resy

    return bbox

Proposed Solution

For now, I will change the get_wgs84_resolution to return a tuple containing width and length of the cell size. But I am still not sure how we decide the cell size of the clipped layer for the case when at least one of them is rectangular.

CC

@timlinux @assefay @Charlotte-Morgan

akbargumbira added a commit to assefay/inasafe that referenced this issue Oct 30, 2014

akbargumbira added a commit to assefay/inasafe that referenced this issue Oct 30, 2014

@timlinux timlinux added the Ready label Nov 4, 2014

@ismailsunni ismailsunni added this to the Version 2.3 milestone Dec 5, 2014

@timlinux timlinux changed the title from Raster Clipping Issue to We should consider cell height and width when clipping rasters Dec 5, 2014

@akbargumbira akbargumbira removed the Ready label Jan 7, 2015

@akbargumbira akbargumbira removed their assignment Nov 16, 2015

@Gustry Gustry modified the milestone: Version 3.0 Sep 1, 2016

@timlinux

This comment has been minimized.

Show comment
Hide comment
@timlinux

timlinux Jan 30, 2017

Contributor

Thankfully another thing we don't need to deal with since we polygonise rasters. @wonder-sk I am closing here but it would be great if you comment on whether the QGIS Raster Alignment tool works well with different X & Y dimensions to raster pixels.

Contributor

timlinux commented Jan 30, 2017

Thankfully another thing we don't need to deal with since we polygonise rasters. @wonder-sk I am closing here but it would be great if you comment on whether the QGIS Raster Alignment tool works well with different X & Y dimensions to raster pixels.

@timlinux timlinux closed this Jan 30, 2017

@wonder-sk

This comment has been minimized.

Show comment
Hide comment
@wonder-sk

wonder-sk Jan 31, 2017

Contributor

Hi there... QGIS raster alignment tool indeed work seamlessly also with non-square cells - both for input rasters and for output ones.

Contributor

wonder-sk commented Jan 31, 2017

Hi there... QGIS raster alignment tool indeed work seamlessly also with non-square cells - both for input rasters and for output ones.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment