[GRASS-user] Converting raster to polygon creates too many polygons

[Markus Metz]

Hanlie Pretoriuswrote:
> Hi,
>
> I've created subcatchments using r.watershed and I would like to
> convert these to polygons. My result at the moment contains too many
> polygons - 194. The original raster had only 9 areas, which should all
> have been contiguous seeing that they came from r.watershed.
>
> Here is the sequence of commands that I use, starting with r.report
> for the information about the subcatchments raster:
> -----
>  |   Type of Map:  raster               Number of Categories: 18
>  |   Data Type:    CELL
>  |   Rows:         33443
>  |   Columns:      33008
>  |   Total Cells:  1103886544
>  |        Projection: Transverse Mercator
>  |            N: -3122229.69982072    S: -3155672.30728643   Res: 0.99998826
>  |            E: -33474.90929218    W: -66483.27873573   Res: 1.00001119
>  |   Range of data:    min = 2  max = 18
> -----

In GRASS, you have full control of the computational region,
consequently you must be able to explain (to the client/supervisor)
every aspect of the computational region. Why is the horizontal
resolution not only slightly off from exactly one meter, but also
different for north-south and east-west? If these are the settings of
the input DEM, why was the DEM created with these settings? If these
resolution settings are different from the DEM's settings, there is
again lots of explaining and justification to do. Life is much easier
with, in this case, exactly 1m resolution which might need to be set
when preprocessing and creating the high-res DEM. The resolution
settings are just slightly different from 1m which can cause rounding
errors later on (export, import, file formats storing resolution in
single and not double precision, other GIS software using single and
not double precision for resolution, etc.)

>
> # Set the resolution to 25m
> g.region rast=c83a_border res=25
Note that this gives you nsres = 24.99447494 and ewres = 25.00634049,
see comment above. Also note that after converting to vector, the
vector will not match the catchment raster because of the changed
resolution. If you are interested in catchments at 25 meter
resolution, resample the DEM to 25 meter, set the computational region
now, don't touch it later, run r.watershed again, convert catchments
to vector.

>
> # Convert the subcatchment raster to vector
> r.to.vect --overwrite
> input=c83a_dem_25m_clipped_basins_100000_clipped at PERMANENT
> output=c83a_subcatchments at PERMANENT feature=area

If you use r.to.vect -v, the categories will correspond to cell
values, making it easier to evaluate the vector.
>
> v.report for the resulting vector gives me 194 categories. The
> original raster had only 9 areas, which should all have been
> contiguous seeing that they come from r.watershed.

The following 5x5 matrix, theoretical example for two basins

1 2 2 2 2
2 1 2 2 2
2 2 1 2 2
2 2 2 1 2
2 2 2 2 2

would be converted not to 2 but to 5 vector areas because the cells
with value 1 are rendered as 4 separate square areas, if you use
r.to.vect -v, with the same category i.e. belonging to the same basin.
BTW, the pattern in the above sample matrix occurs more often with SFD
than with MFD, particularly on such high resolution (roundabout 1m).

Further on, I observed r.to.vect to create redundant areas which can
be removed with v.dissolve.

>
> What am I doing wrong?

Not much apart from messing with the resolution. g.region res=1 -a or
g.region align=c83a_dem_25m_clipped_basins_100000_clipped to cleanly
align a clipped region to a raster map ;-)

> I want an output vector with 9 polygons
> corresponding to the original 9 areas in the subcatchment raster.

You will likely get a vector map with more than 9 areas, but you can
get a vector map with exactly 9 different categories.

Markus M