On 01-12-16 12:01, Moritz Lennert wrote:
On 30/11/16 09:40, Paulo van Breugel wrote:

On Tue, Nov 29, 2016 at 5:44 PM, Moritz Lennert


> An area() function in r.mapcalc would be nice...

Would indeed be a nice addition, but it isn't too difficult to compute,
using:

|g.region rast=afripop10@PERMANENT|
|PI=3.1415926536|
|export `g.region -g`|
|export `g.proj -j`|

This doesn't work for me as the output is

g.proj -j

+proj=longlat
+no_defs
+a=6378137
+rf=298.257223563
+towgs84=0.000,0.000,0.000

and +a is not acceptable as variable name.

You are right, that should  be g.proj -g I think


|r.mapcalc "rcs = (sin(y()+$nsres/2) - sin(y()-$nsres/2)) * ($ewres *
$PI/180) * float($a)^2";|


I tested this approach as well. I.e.

area = above formula
pop_density = pop / area

reproject pop_density

pop = pop_density * (ewres*nsres)

It is faster than the r.in.xyz approach. But it does not seem to be as precise.

I still lost about 100.000 inhabitants, and more when I smooth more (the "loss" increases from nearest neighbor to bilinear to bicubic).

Just guessing, but an issue I had in the past with layers like Afripop is that it may contain cells with extremely high values (a result of the inter/extrapolation methods used to create those layers I guess, not sure it is an issue for e.g., Wordpop). Depending on the resampling method, this could also add to the differences in totals.


In order to avoid precision issues, I tried with both density by m2 and density by km2, but results are similar.
Does it make a difference if you use a much higher target resolution?

I don't know which part of the error comes from the use of spheroid approximation and which part comes from the resampling at reprojection.

It would be interesting to see how results are if using the same routine as above (convert to densities, reproject, convert back to totals), but using the more precise G_area_of_cell_at_row() function to compute the cell areas.


But I do think that if it is important to maintain exact totals, then the r.to.vect - v.proj - v.out.ascii - r.in.xyz solution seems to be more reliable.

Moritz



_______________________________________________
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

Reply via email to