[GRASS-user] how to calculate area of voronoi polygons considering the boundary of the catchment

daljeet sachdev_dj at yahoo.com
Wed Sep 2 15:03:26 EDT 2009


Thanks.

I am using the first method that you suggested. But it fails to create
centroids whereas there are multiple isohyetal strips now and each should
have been given the centroid. To give you the details on what i did.

In real scenario, we will have the defined catchment area (which will be a
close polygon shaped vector) with rainfall measurements on various gauge
points (some gauge points will be within the catchment boundary and some
will be outside that catchment area, but within the region of the catchment
vector). Please correct me if I am wrong.

I am simulating the above scenario by
a) making a xyz point data file. x,y being lat and long and z being the
measurement of rainfall.
b) catchment is represented by a vector (basin_v)

> g.region vect=basin_v

#import the xyz pt data
> v.in.ascii x=2 y=1 cat=0 input="rain_data" out="rainpoint_v" columns="x
> double precision, y double precision, rain_info double precision" fs="\t"
> skip=1 --o

>MASK applied (which represents the raster of the catchment area with each
cell value equal to 1) 

>v.surf.rst rainpoint_v zcol=rain_info elev=rain_rst --o

>r.contour step=1 input=rain_rst output=rain_con_v

At this stage, i have two seperate vectors, a) contour vector (rain_con_v)
which has the lines only b) vector representing the catchment area (basin_v)
whose topo details are as follows

>v.info rain_con_v -t
nodes=7
points=0
lines=7
boundaries=0
centroids=0
areas=0
islands=0
faces=0
kernels=0
primitives=7
map3d=1

Please note the below output. This may be reason of failure.
>v.info basin_v -t
nodes=2
points=0
lines=0
boundaries=1
centroids=1
areas=1
islands=1
faces=0
kernels=0
primitives=2
map3d=0

Now, as you suggested, i need to apply the v.clean (with break) to create
the nodes at the intersection of isolines and boundary of the catchment
vector.

And for the input of v.clean, I patched the above two vectors using v.patch.
Please correct if this is the right step or something else would have been
done?

> v.patch input=rain_con_v,basin_v out=rain_comp_v
WARNING: Vector map <rain_comp_v> already exists and will be overwritten
Patching vector map <rain_con_v at RnD>...
Patching vector map <basin_v at RnD>...
Building topology for vector map <rain_comp_v>...
Registering primitives...
9 primitives registered
3643 vertices registered
Building areas...
 100%
1 areas built
1 isles built
Attaching islands...
 100%
Attaching centroids...
 100%
Number of nodes: 9
Number of primitives: 9
Number of points: 0
Number of lines: 7
Number of boundaries: 1
Number of centroids: 1
Number of areas: 1
Number of isles: 1
Intersections at borders will have to be snapped
Lines common between files will have to be edited
The header information also may have to be edited
v.patch complete. 2 vector maps patched

rain_comp_v looks like this

http://n2.nabble.com/file/n3568825/rain_comp_v.png 

> v.clean in=rain_comp_v out=rain_comp tool=break --o
--------------------------------------------------
Tool: Threshold
Break: 0.000000e+00
--------------------------------------------------
WARNING: Vector map <rain_comp> already exists and will be overwritten
Copying vector lines...
Rebuilding parts of topology...
Building topology for vector map <rain_comp>...
Registering primitives...
9 primitives registered
3643 vertices registered
Number of nodes: 9
Number of primitives: 9
Number of points: 0
Number of lines: 7
Number of boundaries: 1
Number of centroids: 1
Number of areas: -
Number of isles: -
--------------------------------------------------
Tool: Break lines at intersections
--------------------------------------------------
Rebuilding topology for output vector map...
Building topology for vector map <rain_comp>...
Registering primitives...
28 primitives registered
3662 vertices registered
Building areas...
 100%
1 areas built
1 isles built
Attaching islands...
 100%
Attaching centroids...
 100%
Number of nodes: 28
Number of primitives: 28
Number of points: 0
Number of lines: 26
Number of boundaries: 1
Number of centroids: 1
Number of areas: 1
Number of isles: 1

>v.centroids rain_comp out=rain_cen option=add --o
WARNING: Vector map <rain_cen> already exists and will be overwritten
Processing features...
0 new centroids placed in output map
Copying attribute table(s)...
Building topology for vector map <rain_cen>...
Registering primitives...
28 primitives registered
3662 vertices registered
Building areas...
 100%
1 areas built
1 isles built
Attaching islands...
 100%
Attaching centroids...
 100%
Number of nodes: 28
Number of primitives: 28
Number of points: 0
Number of lines: 26
Number of boundaries: 1
Number of centroids: 1
Number of areas: 1
Number of isles: 1
v.category complete. 0 features modified.

Area is still 1 and centroid is 1. If i have understood it correctly, then
we should have multiple area representing the isohyetal strips as output.

Thanks
Daljeet


Micha Silver wrote:
> 
> Hi Daljeet
> 
> 
> daljeet wrote:
> 
>> Hi Micha,
>>
>> I could create isohyetal lines using the suggested steps:
>>
>>   
> ... clipped...
>>
>> Following are the doubts/questions
>>
>> 1. I saw the details of the raster created from v.surf.rst
>>
>>   
>>> r.info rainraster
>>>     
>>
>> zmin_data=0.000000, zmax_data=141.322000                                |
>> zmin_int=-13.588332, zmax_int=52.589905 
>>
>> zmin_data and zmax_data seems to be correct but I do not understand 
>> a) why the interpolated value zmin_int=-13.588332 is -ve
>> b) and what is meaning of zmax_int equal to  52.58990
>>   
> If I understand correctly, the zmax_int and zmin_int are the max and min 
> displacements of the interpolated values vs. the actual data points. The 
> RST method creates an interpolation which does not usually go thru the 
> actual data points, but rather "smooths" over them. These numbers show 
> the drift above and below actual data points. So you have one point 
> where the interpolated rainfall is 13 mm (or whatever the rain units 
> are) below the actual data, and one point where the interpolation is 52 
> mm above.
>> 2. The target is to be able to calculate the area between two isohyetal
>> lines and multiplying that area with the average rainfall within the two
>> isohyets. This is to be done for all isohyetal lines. 
>>
>> Average can be calculated by averaging the value of rainfall of the two
>> isohyetal lines.
>>
>> But how to calculate the area b/w the two isohyetal lines?
>>   
> Well I can think of two approaches to this.
> 1- Using vectors:  You might try to convert the isolines to closed 
> boundaries, then add centroids and find some way to get the rainfall 
> values into the centroids. This would require first to create an "edge" 
> line all around the region (v.in.region will do this) to close the ends 
> of the isolines,  then run v.clean with the "break" tool to get nodes 
> where each contour line crosses this region line. Then you should be 
> able to convert the lines to closed boundaries (v.type) and finally 
> v.centroids to create centroids in each of the resulting strips between 
> isohyetal lines.
> 
> But this procedure looks way to convoluted to me. Instead maybe you 
> could try:
> 2- Using rasters:  You already have the rainfall raster. Perhaps create 
> a reclass raster map of the rainfall raster using discrete values 
> instead of the continuous range. So your reclass rules would look 
> something like:
> 0 thru 20 = 10
> 21 thru 40 = 30
> 41 thru 60 = 50
> ...
> This will group together all values in the interval between 0-20 and 
> give them the value midway between, etc. Then you could do r.to.vect 
> fea=area to create a vector map, and finally v.db.addtable and v.to.db 
> to get the areas of each strip into the vector attribute table. (Note 
> that r.reclass requires an integer raster)
> 
> After that long-winded response, I should add one more comment: Rainfall 
> is not a continuous phenomena across a geographic area. There's no 
> reason to expect that between the 40 mm line and the 60 mm line the 
> average rainfall is 50 mm. If you stay with the original, continuous 
> rainfall raster, and use v.rast.stats then the statistics are totalled 
> on a cell by cell basis, and you avoid the "unnatural" results that 
> might pop up when forcing the data into a framework of arbitrary
> intervals.
> 
> Regards,
> Micha
>>   
> 
> _______________________________________________
> grass-user mailing list
> grass-user at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
> 
> 

-- 
View this message in context: http://n2.nabble.com/how-to-calculate-area-of-voronoi-polygons-considering-the-boundary-of-the-catchment-tp3335923p3568825.html
Sent from the Grass - Users mailing list archive at Nabble.com.


More information about the grass-user mailing list