[GRASS-dev] how to reproject a raster map with absolute numbers without losing data

[Markus Metz]

On Thu, Dec 1, 2016 at 12:01 PM, Moritz Lennert <
mlennert at club.worldonline.be> 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...
[...]
>
>> |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).
>
> In order to avoid precision issues, I tried with both density by m2 and
density by km2, but results are similar.
>
> 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.
>
> 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.

But with this approach you get spatial artefacts. Assume there are 9 cell
in source with

10 10 10
10 10 10
10 10 10

they can become in target e.g.

10  0 20
10  0 20
10 10 10

Total count is maintained, but some target cells might receive no source
cell whereas other target cells might receive more than one source cell.

Markus M
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20161201/1d645ba1/attachment.html>