Opened 6 years ago

Closed 6 years ago

Last modified 5 years ago

#3774 closed defect (fixed)

gdal_rasterize can miss points on the extent edges

Reported by: EliL Owned by: warmerdam
Priority: normal Milestone: 1.9.0
Component: Utilities Version: svn-trunk
Severity: normal Keywords: rasterize, edge, half-pixel-shift
Cc: mattbk, peifer

Description

Running gdal_rasterize -tr 100 100 -a Elevation -l town town.shp townraster.tif converts 12 points to 9 pixels (with values) in the output raster. It should convert it to 11 pixels (with values) in the output raster (1 point has a null value and should not get converted).

town.shp is a point layer with 12 points, 11 of those points have values in -a Elevation (1 is null), and in the output tif, 9 pixels have values. The far south and far east points do not have corresponding pixels with values.

Viewing the results in qgis, the points and pixels do not always line up correctly, perhaps it is the same cause for the extent edges to be missing pixels (i.e. those pixels would be in the next row or column of pixels).

Results should be reproducible with town.shp, attached.

Attachments (1)

town.tar.gz (1.4 KB) - added by EliL 6 years ago.
small point shapefile

Download all attachments as: .zip

Change History (8)

Changed 6 years ago by EliL

Attachment: town.tar.gz added

small point shapefile

comment:1 Changed 6 years ago by mattbk

I can confirm this behavior with the data provided and with my own. Discussion is occurring here. Pixel offset problem can be solved by shifting the raster by 0.5 pixels, but this is not a long-term solution.

comment:2 Changed 6 years ago by mattbk

Cc: mattbk added

comment:3 Changed 6 years ago by peifer

Cc: peifer added

Just to give another example: I would assume that the 5 points in test.csv [1] can be rasterised into a 3 x 3 raster with an extents of [xmin ymin xmax ymax] = [0 0 3 3]. This is obviously not true as one can see at [2]. Shifting the extents of the raster in SE direction by 0.5 pixel seems to be necessary to get the expected result, see [3].

[1]

$ cat test.csv
Longitude,Latitude,Name,Value
0.5,0.5,"First point",1
0.5,2.5,"Second point",2
2.5,2.5,"Third point",3
2.5,0.5,"Fourth point",4
1.5,1.5,"Fifth point",5

$ cat test.vrt
<OGRVRTDataSource>
    <OGRVRTLayer name="test">
        <SrcDataSource>test.csv</SrcDataSource>
        <GeometryType>wkbPoint</GeometryType>
        <LayerSRS>EPSG:4326</LayerSRS>
        <GeometryField encoding="PointFromColumns" x="Longitude" y="Latitude"/>
    </OGRVRTLayer>
</OGRVRTDataSource>


----

[2]

$ gdal_rasterize -l test test.vrt test.tif -a Value -tr 1 1 -te 0 0 3 3 &&
  gdal_translate -of aaigrid test.tif test.asc && cat test.asc
0...10...20...30...40...50...60...70...80...90...100 - done.
Input file size is 3, 3
0...10...20...30...40...50...60...70...80...90...100 - done.
ncols        3
nrows        3
xllcorner    0.000000000000
yllcorner    0.000000000000
cellsize     1.000000000000
 0 0 0
 0 2 0
 0 0 5

----

[3]

$ gdal_rasterize -l test test.vrt test.tif -a Value -tr 1 1 -te .5 -.5 3.5 2.5 &&
  gdal_translate -of aaigrid test.tif test.asc && cat test.asc
0...10...20...30...40...50...60...70...80...90...100 - done.
Input file size is 3, 3
0...10...20...30...40...50...60...70...80...90...100 - done.
ncols        3
nrows        3
xllcorner    0.500000000000
yllcorner    -0.500000000000
cellsize     1.000000000000
 2 0 3
 0 5 0
 1 0 4

comment:4 Changed 6 years ago by peifer

Keywords: half-pixel-shift added
Priority: lownormal
Severity: minornormal

comment:5 Changed 6 years ago by Even Rouault

r22998 /trunk/ (3 files in 3 dirs): gdal_rasterize: fix half pixel shift when rasterizing points; make gdal_rasterize utility increase the computed raster extent by a half-pixel for point layers (#3774)

comment:6 Changed 6 years ago by Even Rouault

Milestone: 1.9.0
Resolution: fixed
Status: newclosed

comment:7 Changed 5 years ago by warmerdam

I have backported the llrasterize.cpp portion of this to 1.8 in time for 1.8.2.

Note: See TracTickets for help on using tickets.