On Thu, Nov 19, 2009 at 04:30:55PM +0100, Achim Kisseler wrote: > Hi, > > > use r.to.vect to map every cell to a polygon, project it to the lambert > > map, do the calculation on the polygon, reproject it back to the > > just run v.to.db -p option=area map=... > > then you get the size for every lines cell. > > To get the raster-map: > r.mapcalc "ns_col=if(col()==1 , row() , null())"
My idea was to do the other way around. In the lonlat location (cru_maps), I do a polygon that spans the whole region and delimitates cells: rows_nr=`g.region -g |grep '^rows=' |sed 's/^rows=//'` r.mapcalc "grid_map=col() + ((row()-1) * $rows_nr)" r.to.vect input=grid_map out=grid_map_v feature=area Then in the Lambert map, I project and use v.rast.stats: v.proj input="grid_map_v" location="cru_maps" v.rast.stats vector=grid_map_v raster=lambert_map colprefix=stat then, in the lonlat location I project back v.proj input="grid_map_v" location="network_af" output=result_v v.to.rast input=result_v output=result_rast column=stat_sum Right now I am stuck at the first projection since it triggered an assertion... -- Pat _______________________________________________ grass-user mailing list [email protected] http://lists.osgeo.org/mailman/listinfo/grass-user
