[GRASSLIST:397] Re: a huge matrix ?

orkun temiz at deprem.gov.tr
Mon Jun 16 04:24:54 EDT 2003


>I'm not sure that you found an answer to this yet, Ahmet - if you did, 
>please ignore this suggestion:
>
>1) use r.mapcalc to create an integer combined index for the combinations 
>of your raster layers:
>
>    tmp1 = mapa*10 + mapb
>    tmp2 = tmp1*10 + mapc
>    ....
>
>or some such way of generating an integer index differentiating the 
>different types. If there are raster values outside 0-9, maybe reclass 
>first, or use 100 instead of 10 to left shift the index where there are 
>10-99 discrete values.
>
>2) In R, generate the same integer index, and make a data frame with the 
>integer index in one column, and the probabilities in another.
>
>3) Use the R/GRASS interface (or r.stats -1 and scan()) to move either the 
>whole of your index raster, or more likely tiles made by changing the 
>GRASS region, to R.
>
>4) in R, use vectorised lookup to search for the probability value given 
>the index value - match() is fast and powerful (I haven't checked this, 
>but do use match() in other code - it will probably work much better on 
>tiles than your big raster) and assign the probability values to a new 
>vector.
>
>5) move this vector back to GRASS as FCELL or DCELL. If you are tiling, 
>finish off by patching the tiles back together.
>
>6) My guess is that the tiling, import, matching, and export can all be 
>done in a loop within R.
>
>Rationale: the combined index avoids all the if statements - 11K if's is 
>rather a lot! I did think about pushing it all out say to postgresql and 
>using SQL (one table of RHS raster values, one table of lookups) but there 
>were still seven fields on each side to index on, so it seemed at least 
>worth trying to combine them automatically by left shifting (probably not 
>the right word) and addition, which are cheap operations in r.mapcalc.
>
>Please let me know how you tackled this problem, your solution will have 
>general interest, and should be in the general GRASS knowledge base. I can 
>expand on this off-list if it is of any use.
>
>Roger
>
>
>On Wed, 11 Jun 2003, Glynn Clements wrote:
>
>  
>
>>orkun wrote:
>>
>>    
>>
>>>>>I have a huge matrix containing 7 column and 10500 rows.
>>>>>Each columns stand for different maps. Each rows represent "unique"
>>>>>combinations of catagories of maps and probability values
>>>>>of each row.
>>>>>Something like:
>>>>>
>>>>>             map1  map2 map3 ... prob
>>>>>1   cat1       x     x    x        x
>>>>>2   cat2       x     x    x        x
>>>>>.     .
>>>>>.     .
>>>>>10500 .
>>>>>
>>>>>So I have 10500 if-clause to create a probability map.
>>>>>I tried r.mapcalc (as from a text file) but grass gave "stack overflow" 
>>>>>error.
>>>>>   
>>>>>
>>>>>          
>>>>>
>>>>I'm not surprised.
>>>>
>>>> 
>>>>
>>>>        
>>>>
>>>>>Is there anyone making suggestions ?
>>>>>   
>>>>>
>>>>>          
>>>>>
>>>>You need to describe your problem in more detail; it isn't clear what
>>>>you are trying to calculate.
>>>>
>>>> 
>>>>
>>>>        
>>>>
>>>Thank you for your interest
>>>
>>>You wanted to know the problem in more detail:
>>>
>>>what I did:
>>>1- using r.stat command I produced a text file for statistical calculations
>>>2- using R statistical package I produced probability results of each 
>>>map's for each categories.
>>>it means 7 columns of maps 10500 rows plus probability of each 
>>>combination.   Actually it is a contingency table. Any data in table 
>>>(matrix) is resulting  from unique combination of maps and their
>>>categories.
>>>3- I want to port these results to grass to create a  map (of landslide 
>>>susceptibility map) . Firstly I thought I used r.mapcalc with piping 
>>>from text file. But because text file contains a huge matrix that is 
>>>made up of 10500 if statement (for same  number of  rows) process failed.
>>>I attached a part of text file after modifying for r.mapcalc. you  see 
>>>all ifs are unique.
>>>      
>>>
>>OK, so: each line is a rule; for each pixel you scan the list of rules
>>to find the first match, then output the associated result. Right?
>>
>>The most straightforward approach (although not particularly
>>efficient), would be to use an iterative approach, i.e.
>>
>>	n=0
>>	last='null()'
>>	while read a1 a2 a3 a4 a5 a6 a7 a8 ; do
>>		result=pr2.pass$n
>>		r.mapcalc "$result = if (f28geo5==$a1 && f28slpcat==$a2 && f28concavity==$a3 && n_facing==$a4 && road_bas==$a5 &&str_buf==$a6 && f28f4==$a7, $a8, $last)"
>>		last=$result
>>	done < rules.txt
>>
>>where rules.txt looks like:
>>
>>	21 6 1 1 1 1 0 701.22
>>	 7 4 2 1 1 2 1 702.06
>>	...
>>
>>Alternatively, you could combine this with your initial approach, only
>>matching a limited number of rules on each pass, and using several
>>passes.
>>
>>However, I note that there seem to be very few possible values for
>>most of the maps. In that case, a far more efficient approach would be
>>to:
>>
>>1. Reclass each input map so that all of the significant categories
>>are consecutively numbered from zero (0,1,2,3,...); any insignificant
>>categories would be merged into the last category.
>>
>>2. Compute a composite map, using e.g.
>>
>>	composite = f28geo5 + 22 * (f28slpcat + 9 * (f28concavity + 3 * ...
>>
>>where the numbers are the total number of categories in the preceding
>>map.
>>
>>3. Assign the values as labels for the composite categories.
>>
>>4. Use "@map" in r.mapcalc to use the label instead of the category.
>>
>>
>>    
>>
>
>  
>
Thank you sir

You know my study from previous messages.
I had several good suggestions from so many kind people including you.
I am going to try them especially grass/R interface.
Firstly I want to know if it is correct to use the formulas in the 
attachment.
I divided the file to 7 according to probability values to handle with 
r.mapcalc.
(Beforehand I had sorted the data according to probability vales by a 
spreadsheet). Then I used
r.mapcalc as in attachment. [I  would reclassify in the grass anyway 
according to
probability values as high, moderate and low . Because of a large text 
file . I did division
(reclassifying) in input level ( or in text file level)].
 But I hesitated on addition of each if expression in r.mapcalc. It must 
be okey. Because
each if expressionsare unique, they represent different areas. So added 
values don't coincide each other.

Could you please tell me firstly whether my approach is correct or 
false. And why ?


kind regards


______________________________________
Scanned and protected by Inflex
http://pldaniels.com/inflex

______________________________________
The views and opinions expressed in this e-mail message are the sender's own
and do not necessarily represent the views and the opinions of Earthquake Research Dept.
of General Directorate of Disaster Affairs.

Bu e-postadaki fikir ve gorusler gonderenin sahsina ait olup, yasal olarak T.C.
B.I.B. Afet Isleri Gn.Mud. Deprem Arastirma Dairesi'ni baglayici nitelikte degildir.

-------------- next part --------------
A non-text attachment was scrubbed...
Name: pr2.gz
Type: application/x-gzip
Size: 925 bytes
Desc: not available
Url : http://lists.osgeo.org/pipermail/grass-user/attachments/20030616/c43419e5/pr2.gz


More information about the grass-user mailing list