[GRASS-user] SOLVED: How to implement 2D lookup table?

[Rainer M. Krug]

Moritz Lennert <mlennert at club.worldonline.be> writes:

> On 23/05/13 11:04, Rainer M. Krug wrote:
>>
>> Moritz Lennert<mlennert at club.worldonline.be>  writes:
>>
>>> On 23/05/13 10:26, Rainer M Krug wrote:
>>>>
>>>> Hi
>>>>
>>>> I have two raster layers
>>>>
>>>> ,----
>>>> | dc: CELL, 1, ..., 7, NULL
>>>> | ageClass: CELL, 1,2,3,NULL
>>>> `----
>>>>
>>>> And I want to create a third one, depending on the value in a lookup
>>>> table. dc corresponds to the row, ageClass to the column,
>>>>
>>>> The lookup table looks as follow:
>>>>
>>>> ,----
>>>> |     seedling young adult
>>>> | dc1   425000 87500 37500
>>>> | dc2   270000 55000 23500
>>>> | dc3   127500 26000 11250
>>>> | dc4    32500  7100  3200
>>>> | dc5     5550  1220   500
>>>> | dc6      605   132    55
>>>> | dc7       56    12     6
>>>> `----
>>>>
>>>> Now I want to have something like (using r.mapcalc as an example)
>>>>
>>>> ,----
>>>> | r.mapcalc expression = "result = lookupTable[col=ageClass, row=dc]"
>>>> `----
>>>>
>>>> and e.g. for one cell
>>>>
>>>> ,----
>>>> | ageClass = 3
>>>> | dc = 5
>>>> | ==>
>>>> | result = 500
>>>> `----
>>>>
>>>> I haven't found anything like this in GRASS.
>>>>
>>>> I could use MASKS and then r.reclass(), or generate dynamically a
>>>> complex if() construct, but both options sound a little bit
>>>> awkward to me.
>>>>
>>>> I am scripting GRASS from R, so there is always the option of doing it
>>>> in R, but there is the issue of reading the data into R and writing it
>>>> back to GRASS, which takes time. As I am doing this in a simulation
>>>> which runs several times, this is not an unimportant issue.
>>>>
>>>> Is there better way of achieving this?
>>>
>>> I don't see a possibility to include the lookup table into r.mapcalc,
>>
>> Me neither - that was just the easies to show what I want to do.
>>
>>> but anyway, IIUC, the result you are looking for is not a raster map,
>>> but a numerical value, or ?
>>
>> No - the lookup would be for each cell, therefore resulting in a CELL
>> raster map.
>>
>>>
>>> My best first guess is:
>>>
>>> r.stats -c in=dc,ageClass
>>>
>>> possibly combined with a little awk to do the lookup:
>>>
>>> r.stats -c in=dc,ageClass | awk '{if($1==5&&  $2==3) print $3}'
>>>
>>> Or, if the maps are heavy, feed the data resulting from the r.stats
>>> call into a new db table and then use sql queries to do the lookup. Or
>>> feed the results of r.stats into an R dataframe.
>>
>> Hmm - sounds interesting, but also for other purposes. But as the result
>> will be a raster map, I don't think this approach is usable here (unless
>> I can feed the numbers back into a raster map, i.e. the inverse of
>> r.stats -c)
>
> r.cross + r.reclass ?

There is one aspect in r.cross which makes this approach quite
difficult to implement: the value of the result map layer can not easily
be traced back to the values of the two input layers:

,----
| Note: The actual category value assigned to a particular combination in
| the result map layer is dependent on the order in which the combinations
| occur in the input map layer data and can be considered essentially
| random. The example given here is illustrative only.
`----

But this gave ma an idea: As I know the *theoretical* maximum values of
each map (as defined in the lookup table), I can do the following:

,----
| r.mapcalc expression="lookup = dc + (ageClass-1)*7"
`----

and flatten in the same way the lookup table 
and convert it to reclass
rules:

,----
| 1  =  425000 dc1 seedling
| 2  =  270000 dc2 seedling
| 3  =  127500 dc3 seedling
| 4  =   32500 dc4 seedling
| 5  =    5550 dc5 seedling
| 6  =     605 dc6 seedling
| 7  =      56 dc7 seedling
| 8  =   87500 dc1 young
| 9  =   55000 dc2 young 
| 10 =   26000 dc3 young
| 11 =    7100 dc4 young
| 12 =    1220 dc5 young 
| 13 =     132 dc6 young
| 14 =      12 dc7 young
| 15 =   37500 dc1 adult
| 16 =   23500 dc2 adult
| 17 =   11250 dc3 adult
| 18 =    3200 dc4 adult
| 19 =     500 dc5 adult
| 20 =      55 dc6 adult
| 21 =       6 dc7 adult
`----

Now I can use 

,----
| r.reclass input=lookup output=finalLookup rules=TheRulesAsAbove
| r.mapcalc "finalLookup = finalLookup"
| r.remove rast=lookup
`----

and then I should have it.

Thanks a lot,

Rainer

>
> moritz
> _______________________________________________
> grass-user mailing list
> grass-user at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
>


-- 
Rainer M. Krug
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 489 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20130523/7aa38a0b/attachment-0001.pgp>