Nikos A:
(MM & ML, I really appreciate the time you spend on this)
[...]
AFAIU, MarkusM's change only affects dmax>0, so it's normal that you
don't see a difference.
Yep, it's clear.
> [2] map #2 v.info -t pareto_...<Very Long Name>
> nodes=26956 points=0 lines=0 boundaries=26955 centroids=6650
> areas=6650 islands=1 primitives=33605 map3d=0
So, you have 6000+ areas times 26000+ points.
Well, this is just a "spearfish" small example. Worst case beforehand here is:
areas=2376 _vs_ points=593273 (in a region of rows:1201, cols:550,
cells:660550 -- but this does not make any difference though)
> Hmmm... ? So, if use my 600000+ points maps, you can imagine why I
> report timings>20h.
You mean that it's normal that it takes so long, or that it is not
normal, and that's why we can imagine that you report it ? 
It was a (very) wild guess, smthn like "3m5.339s for 23.000, how much for 600K
thinking that time-consumation does _not_ increase linearly, but more
aggressive).
Here's an attempt with a 100x100 grid (=10000 areas) and the 600000
points randomly created:
v.mkgrid mygrid n=100,100 --o
Can someone confirm that even this takes (a bit) long (or is it only my box)?
time v.distance from=mypoints@sqlite to=mygrid upload=to_attr
column=cat_municip to_column=cat dmax=0.0
real 2m22.505s
GRASS 6.4.0svn (nc_spm_08):~ > real 1m42.860s [ 
| 
]
(shrug)
Don't know what's wrong here!? To make the long story short following my
scripts:
1. pareto_1 converts all pixels of Landsat7 samples to (centroids ->) points.
Yes, res=30m, many pixels (e.g. the 593273 above). Includes only: g.region,
r.to.vect, v.db.dropcol.
2. pareto_2 loops over the point(s) maps and creates "low-resolution" vector
grid (using v.mkgrid, e.g. 2376 boxes that correspond to a raster with
res=500m). Includes only: g.region, v.in.region, v.select, g.remove, v.mkgrid
3. pareto_3 just counts points in polygons. Includes: v.db.addcol, v.what.vect
(this is the time-killer), v.support, v.db.addcol, v.db.update (also
inefficiently executed but not the real problem)
(...the rest of pareto_? is another story...)
I will give another go with the spearfish60 dataset but this time keeping the
high-resolution to 30m and the low resolution to 500m. In detail, all of the
commands that have to be executed are:
--%<---
grass64 /spearfish60/user1/
g.copy rast=landcover.30m,landcover.30m
g.region rast=landcover.30m -pa # cells: 294978
r.mapcalc "rangeland.30m = if((landcover.30m == 51 || landcover.30m == 71 ||
landcover.30m == 81 || landcover.30m == 92), 2, null())"
r.grow input=rangeland.30m output=pareto_ref_30m radius=1.1 metric=manhattan
r.mapcalc "coi_1 = if((landcover.30m == 51 || landcover.30m == 71 ||
landcover.30m == 81 ), 2, null())"
r.mapcalc "coi_2 = if(( landcover.30m == 71 || landcover.30m == 92), 2,
null())"
g.region rast=landcover.30m res=500 -pa
r.grow input=coi_1 output=pareto_classification_500m_1 radius=1.1
metric=manhattan
r.grow input=coi_2 output=pareto_classification_500m_2 radius=1.1
metric=manhattan
g.mremove rast=coi_[12] -f
pareto_1_vectorise_rasters.py reference_raster=landcover.30m
reference_coi_rasters=pareto_ref_30m
classification_rasters=pareto_classification_500m_1,pareto_classification_500m_2
pareto_2_create_lowres_vector_grid.py highres=30 lowres=500 --v
# here would go the "pareto_3" script, running many times v.what.vect
depending on the number of "area-maps" and naturally the "poin-maps". However,
one test-command is enough to demonstrate:
time v.what.vect --v vector=pareto_ref_points___pareto_ref_30m
qvect=pareto_grid___OVER___ref_coi_points___pareto_grid___OVER___ref_30m___ON___classification_V___classification_1_500m
column=gridcell_1 qcolumn=cat
-->%---
Now this is crawling and I am waiting [*]. Didn't wanna let you waiting too
long to know that this keeps me waiting - so I've hit the Send button :-p
Nikos
---
[*] You Waited Too Long, Casey Bill Weldon (1936)
